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ABSTRACT 

We present an N-body model that reproduces the morphology and kinematics of the Magellanic 
Stream (MS), a vast neutral hydrogen (HI) structure that trails behind the Large and Small Magellanic 
Clouds (LMC and SMC, respectively) in their orbit about the Milky Way. After investigating 8 x 10 6 
possible orbits consistent with the latest proper motions, we adopt an orbital history in which the 
LMC and SMC have only recently become a strongly interacting binary pair. We find that their 
first close encounter ^2 Gyr ago provides the necessary tidal forces to disrupt the disk of the SMC 
and thereby create the MS. The model also reproduces the on-sky bifurcation of the two filaments of 
the MS, and we suggest that a bound association with the Milky Way is required to reproduce the 
bifurcation. Additional HI structures are created during the tidal evolution of the SMC disk, including 
the Magellanic Bridge, the "Counter-Bridge," and two branches of leading material. Insights into the 
chemical evolution of the LMC are also provided, as a substantial fraction of the material stripped 
away from the SMC is engulfed by the LMC. Lastly we compare three different N-body realizations 
of the stellar component of the SMC, which we model as a pressure-supported spheroid motivated 
by recent kinematical observations. We find that an extended spheroid is better able to explain the 
stellar periphery of the SMC, and the tidal evolution of the spheroid may imply the existence of a 
stellar stream akin to the gaseous MS. 

Subject headings: Magellanic Clouds — Galaxy: halo — galaxies: kinematics and dynamics — galaxies: 
evolution 



1. INTRODUCTION 

Disentangling the interaction history of the Large and 
Small Magellanic Clouds (LMC and SMC, respectively; 
"MCs" collectively) has been facilitated by modeling the 
formation of the Magellanic Stream (MS), a prominent 
neutral hydrogen (HI) structure that trails behind the 
MCs in their orbit about the Milky Way (MW). Theoret- 
ical models have assigned one of two physical mechanisms 
for creating the MS: tidal stripping (Lin & Lynden-Bell 
1977; Murai & Fujimoto 1980, MF80; Gardiner, Sawa & 
Fujimoto 1994, GSF94; Gardiner & Noguchi 1996, GN96; 
Yoshizawa & Noguchi 2003; Connors et al. 2006, C06) 
and ram pressure stripping (Meurer et al. 1985; Heller 
& Rohlfs 1994; Moore & Davis 1994; Mastropietro et al. 
2005, M05). The common feature of these MS formation 
models is their reliance on multiple strong interactions 
between the MCs and the Milky Way. The plausibil- 
ity of such an interaction history for the MCs has been 
challenged, however, by the high-precision proper mo- 
tion measurements of Kallivayalil et al. (2006a,b; K06). 
These proper motion estimates imply such large orbital 
velocities for the MCs that a "first passage" scenario has 
been proposed in which the MCs are not bound to the 
Milky Way but are instead passing by for the first time 
(Besla et al. 2007). Because the past orbits of the LMC 
and SMC have been called into question, the origin of the 
MS must also be revisited. For example, the "blowout 
hypothesis" has emerged recently as an alternate sce- 
nario for the formation of the MS, which relies on inter- 
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nal mechanisms to expel gas from the MCs (Olano 2004, 
Nidever et al. 2008). 

Attempts to reconcile the K06 proper motions with the 
formation of the MS have not yet produced a compelling 
dynamical model. For example, even though the first 
passage model of Besla et al. (2010, B10) reproduces a 
tidal tail resembling the MS, the model does not ade- 
quately recover the known positions and velocities of the 
MCs. In particular, a number of important orbital con- 
straints for the SMC are not satisfied, including on-sky 
positio n, li ne-of-sight velocity, and proper motion (see 
section |B~4| for a clarifying discussion as well as Besla et 
al. 2012). Ruzicka et al. (2010) searched a large param- 
eter space centered on the K06 proper motions, but they 
were unable to find a test-particle model that could si- 
multaneously explain the on-sky location and kinematics 
of the MS. In our previous work, we were able to identify 
a good MS model, but we had to adopt a massive isother- 
mal halo for the Milky Way in order to combat the large 
velocities of the MCs (Diaz & Bekki 2011a, DBlla). 

Nevertheless, one may salvage some insight from the 
above three studies. Ruzicka et al. (2010) point out 
that their most promising models exhibit a common in- 
teraction history: two close encounters between the LMC 
and SMC at <2.5 Gyr and 150 Myr ago, each of which 
triggers an epoch of tidal stripping from the SMC disk. 
Similar conclusions are reached in DBlla and B10, in 
which the MS forms during the first strong binary inter- 
action between the MCs rather than during an interac- 
tion with the Milky Way. Whereas previous tidal mod- 
els relied on a combination of tides from both the Milky 
Way and LMC (e.g., GN96, C06), these recent models 
indicate that the onset of violent LMC tidal interactions 
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could have been sufficient to disrupt the SMC disk and 
thereby create the MS. This new scenario would imply 
that the binary state of the MCs is a recent phenomenon 
and that the MS is an artifact from their dynamical cou- 
pling. Furthermore, observations appear to corroborate 
this scenario, particularly because the star formation his- 
tories of the LMC and SMC exhibit two correlated bursts 
of star formation at ~ 2 Gyr ago and ~500 Myr ago (Har- 
ris & Zaritsky 2009). These two epochs may correspond 
to two strong tidal interactions suffered mutually by the 
LMC and SMC, which in turn would create the physical 
conditions necessary for episodic bursts of cluster forma- 
tion (Bekki et al. 2004; Piatti et al. 2005; Piatti 2011). 

Regarding the tension between the K06 proper motions 
and MS formation models, a possible compromise has 
been offered by the recent proper motion measurements 
of Vieira et al. (2010, V10) and Costa et al. (2009): 
perhaps the velocities of the MCs are indeed larger than 
traditionally assumed (e.g., C06, M05) but not as large 
as what K06 would imply. The V10 measurements are 
not nearly as precise as those of K06, but Bekki (2011a) 
has suggested that the V10 proper motions may never- 
theless be more accurate than those of K06 due to the 
larger sample size (3822 stars versus 810 stars, respec- 
tively, for the LMC) and the longer baseline (40 years 
versus 2 years, respectively). Bekki (2011a) argues that 
random stellar motions would create unknown system- 
atic errors that could undermine the accuracy of the K06 
measurement but would be largely suppressed for V10. 
Another unique advantage of V10 is that they measure 
a precise relative proper motion between the LMC and 
SMC, achieved by tracking the MCs in the same wide 
field images (~450 square degrees). Considering that 
many properties of the MCs can be understood in terms 
of their recent activity as a binary pair (star formation 
histories, Harris & Zaritsky 2009; bursts of star cluster 
formation, Bekki et al. 2004, Piatti et al. 2005; for- 
mation of MS, DBlla, B10), the V10 constraint on the 
relative motions of the LMC and SMC are an asset to 
theoretical models. 

In the present work we use orbital models and N- 
body simulations to address the formation of the MS 
and other HI structures of the Magellanic system, in- 
cluding the Magellanic Bridge, which extends between 
the SMC and LMC, and the Leading Arm (LA), which 
stretches ahead of the MC orbits. The framework of 
the present model provides an improvement over previ- 
ous MS models (e.g., C06, M05, B10) in three important 
ways. First, we explore a wide range of orbital histo- 
ries with velocities constrained by the observational re- 
sults of V10. Second, we adopt a multi-component po- 
tential for the Milky Way composed of a disk, bulge, 
and NFW halo, which is a more realistic choice than, 
for instance, the one-component Milky Way potentials 
adopted by C06, B10, and DBlla. And third, we rep- 
resent the SMC as having both a rotating disk and a 
non-rotating spheroid, whereas previous MS models have 
only represented the SMC as a "pure disk" system. We 
are thus able to correctly reproduce the observed kine- 
matics of the SMC, both its rotating HI component (Sta- 
nimirovic et al. 2004) and its non-rotating stellar com- 
ponent (Harris & Zaritsky 2006, HZ06). Our adopted 
spheroid model is also consistent with the extended stel- 
lar halo of the SMC (e.g., De Propris et al. 2010; Nidever 



et al. 2011). 

In our best model, we find that the SMC disk remains 
intact until the dynamical coupling of the MCs ~2 Gyr 
ago, at which point the MS and LA are violently torn 
away by the strong tidal forces of the LMC. Additional 
tidal debris is subsequently engulfed by the LMC, creat- 
ing a transfer of mass from the SMC that may provide in- 
sight into the chemical enrichment history of the LMC. A 
second tidal encounter occurring ^250 Myr ago is respon- 
sible for pulling the Bridge from the SMC disk as well 
as a complementary structure that we call the "Counter- 
Bridge" . Our model exhibits strong agreement with the 
HI kinematics and morphology of the MS, particularly its 
well-known bifurcation into two distinct filaments (Put- 
man et al. 2003a). For the first time, we are able to 
provide a structural interpretation of the MS filaments 
within the context of a dynamical model (however, see 
Nidever et al. 2008), and we find that a bound associ- 
ation with the Milky Way is required for at least 2 Gyr 
to explain the bifurcation. We accordingly suggest that 
the morphology of the MS provides a strong argument in 
favor of bound orbits for the MCs, and in particular that 
a first passage orbit, in which the MCs interact with the 
Milky Way only recently, may not be able to explain the 
MS bifurcation. 

Using "disk plus spheroid" N-body models, we describe 
for the first time the evolution of the SMC spheroid un- 
der the same tidal forces that pulled the MS from the 
SMC disk. In order to provide a comparative anal- 
ysis of important physical parameters, three different 
spheroid models are presented. We demonstrate that an 
extended, low-density spheroid is better able to repro- 
duce observations, including the stellar kinematics of the 
SMC (HZ06), the recent discovery of a break population 
of red giants (Nidever et al. 2011), and the observation 
of kinematically and chemically peculiar stars within the 
LMC (Olsen et al. 2011). As a consequence of reproduc- 
ing these observations, we find that the tidal evolution of 
the extended spheroid necessarily predicts that a coher- 
ent stellar stream is stripped away under the same forces 
that remove the MS and LA from the disk. The on-sky 
location of the predicted stellar stream is slightly offset 
from the observed MS, particularly at its tip, which pos- 
sibly explains why previous attempts to observe a stellar 
counterpart to the MS have failed (e.g., Guhathakurta & 
Reitzel 1998). 

The plan of the paper is as follows: in the next section, 
we outline our numerical model including a brief discus- 
sion of adopted parameters. We represent the SMC as a 
multicomponent system (disk, spheroid, and dark mat- 
ter halo), and we describe the evolution of its disk and 
its spheroid in separate sections for the sake of clarity. 
In section [3l we present our results concerning the tidal 
evolution of the SMC disk, including a discussion of the 
MS bifurcation. In section H] we compare three different 
spheroid models which enables us to discuss the possibil- 
ity of a tidal stream of stars. In section [5] we provide a 
general discussion including a comparison of recent MS 
formation models. In section |6] we summarize and con- 
clude. 

2. THE NUMERICAL MODEL 

The present investigation is three-fold. First, we search 
for the most plausible and realistic orbits of the MCs with 
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Fig. 1. — The proper motion of the LMC (left) and SMC (right). Observational ellipses outline the 68.3% confidence regions and la error 
bars for Vieira et al. (2010; V10, solid) and Kallivayalil et al. (2006a, b; K06, dashed). The parameter space of possible values explored 
in the present work are shaded by the grey boxes, and the adopted values are indicated by the star. The values chosen by Gardiner & 
Noguchi (1996; plus sign) and Diaz & Bekki (2011a; cross) are also shown for comparison. Models that fail to reproduce the SMC orbit, 
e.g., Mastropietro et al. 2005 and Besla et al. 2010, are not included in the plot (see text). The Vieira et. al (2010) proper motion for 
the SMC is measured relative to the LMC, and the corresponding ellipse therefore roams over the large grey box (right panel) as the LMC 
proper motion is varied. For the purposes of the figure, we have chosen to fix the ellipse by the LMC proper motion adopted in the present 
model (star, left panel). 



respect to the Milky Way by using a backward integra- 
tion scheme (MF80; GSF94; GN96). We run - 8 x 10 6 
orbital models based on the full range of our parame- 
ter space, and we then focus our attention on the subset 
of models in which the MCs form a recent binary pair. 
The recent formation of a strong binary state is assumed 
to be necessary for the formation of the MS (DBlla). 
Second, we run a host of test particle simulations and 
compare the spatial distributions against that of the ob- 
served MS. This allows us to further narrow the set of 
promising orbits. And third, we investigate the struc- 
ture and kinematics of the simulated MS using medium 
resolution N-body models. In this last phase of the in- 
vestigation, we are able to isolate the best orbital model 
within our parameter space. We present our findings as 
a set of high resolution N-body simulations that explain 
the observed properties of the MS, LA, Bridge, and stel- 
lar structure/kinematics of the SMC in a self-consistent 
manner. 

2.1. Orbital Models 

Despite the fact that many proper motion studies have 
been carried out in the last two decades (e.g., see V10 
for discussion), an accepted set of measured values has 
not been converged upon. Accordingly, we constrain our 
orbital models by the on-sky positions, distances, and 
radial velocities of the MCs (see Table fl] for adopted val- 
ues) whereas we explore a large possible range of proper 
motions. We center our investigation on the measured 
values of V10 for the absolute proper motion of the 
LMC, (^ Q cos<5, w) lmc = (1.89 ±0.27, 0.39 ±0.27) mas 
yr _1 and for the more precise relative proper motion of 
the SMC about the LMC (/i Q cos<S, ^)smc = (^ Q cos5, 
M<5)lmc - (0.91 ± 0.16, 1.49 ± 0.15) mas yr" 1 . Figure p] 
gives the investigated range as grey boxes for the LMG 



(left panel) and for the SMC (right panel). Notice that 
as we select different values for the LMC proper motion, 
the corresponding V10 constraint on the SMC (solid el- 
lipse, right panel) roams across the plane to different 
values. The range of proper motions that we consider 
is sufficiently large to overlap with the la errors of the 
K06 measurement, shown as a dotted ellipse for the LMC 
(left) and SMC (right). To avoid clutter in the figure we 
plot only the observed values of V10 and K06, but other 
recent measurements deserve to be mentioned such as 
Piatek et al. (2008) and Costa et al. (2009). 

The adopted proper motion values of the LMC and 
SMC in the present work are indicated by the stars in 
Figure [T] We selected these particular values after con- 
sidering ~ 8 x 10 6 orbital models for the MCs (see below 
for details). The adopted values fall between the pre- 
viously assumed values of GN96 (plus sign) and DB11 
(cross). In this sense, the present work is a compro- 
mise between the classic bound scenario for the MCs 
(i.e., GN96) and the high velocity orbits implied by K06. 
There are a host of tidal models for the formation of the 
Magellanic Stream which use the same proper motions 
assumed by GN96, including GSF94, Gardiner (1999), 
Yoshizawa & Noguchi (2003), and Connors et al. (2006). 
The plus sign in Figure [l] represents all of these previ- 
ous works. Models that reproduce only the LMC orbit 
have been excluded from the plot, including the ram pres- 
sure stripping model of Mastropietro et. al (2005), which 
omits the SMC entirely, and the first passage model of 
Besla et al. (2010), which does not reproduce the correct 
SMC position or velocity (see Besla et al. 2012). 

We calculate the past orbits of the LMC and SMC with 
respect to the Milky Way by adopting the backward in- 
tegration scheme originally devised by MF80. To carry 
out these orbit integrations, we must assume model pa- 
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TABLE 1 

Mass and Orbital Parameters 



Property 



Adopted Value 



Milky Way 



t Total Mass (r <300 kpc) 1 ' 2 
Disk Mass (Mj) 1 
Bulge Mass (Mb) 1 
NFW Virial Mass (M vir ) 2 
NFW Virial Radius (R vil ) 
Circular Velocity (V c i r ) 3 
*Distance to Sun (i? ) 3 ' 4 
^Velocity of Sun 5 



1.73 x 10 12 M * 
5.0 x 10 10 M 
0.5 x 10 10 M 
1.30 x 10 12 Mq* 
175 kpc* 
240 kms' 1 * 
8.5 kpc 

(10.0, 5.2 + V cir , 7.2) km s" 1 



LMC 



SMC 



Mass (M mc ) b ' 7 
Scale radius (a mc ) 
Proper Motion (/i a cos S, /is) 8 
Line-of-Sight Velocity 9 ' 10 
On-Sky Position (a, <5) 9,11 
Distance Modulus 12 ' 13 
tPosition Vector {X ,Y ,Z) 
*Space Velocities (U,V,W) 



(1 



10 10 M 
3 kpc 
, 0.30) mas yr — 1 
262 km s" 1 
(81.9°, -69.9°) 
18.50 

(-0.8, -41.6, -27.0) kpc 
(-50.7, -226.1, 229.3) km s" 



3 X 10 9 Mq 

2 kpc 

(0.92, -1.37) mas yr" 1 
146 km s" 1 
(13.2°, -72.5°) 
18.95 

(15.3 -36.9 -43.3) kpc 
(-4.2, -223.5, 191.0) km s" 



1 Binney & Tremaine (2008); 2 Gnedin et al. (2010); 3 Reid et al. (2009); 4 Gillessen et al. (2009); 
5 Dehnen & Binney (1998); 6 Kim et al. (1998); 7 Stanimirovic et al. (2004); 8 Vieira et al. (2010); 
9 van der Marel et al. (2002); 10 Harris & Zaritsky (2006); "Piatek et al. (2008); 12 Freedman et 
al. (2001); 13 Cioni et al. (2000) 

t Includes bulge mass, disk mass, and NFW halo mass within r < 300 kpc. 
t Measured in a galactocentric frame of rest. 

* We find that the tidal evolution of the SMC is nearly identical to the present model when the 
following parameters are substituted for the NFW halo: M vir =1.90 X 10 12 M , R vir =269 kpc. 
In this "alternate model" (see text), the total mass of the Milky Way within 300 kpc is 2.06 X 
10 12 M and the circular velocity at the solar radius is 220 km s — 1 . The parameters for the 
LMC and SMC are same as above, expect for the V component of the space velocities which 
decrease by 20 km s" 1 owing to the change in circular velocity. 



rameters for the following quantities: (1) the shape of 
the Milky Way's gravitational potential as a function of 
distance r from the Galactic center, (2) gravitational po- 
tential of the MCs, (3) total masses of the MCs, and (4) 
the form of dynamical friction between the Milky Way's 
dark matter halo and the LMC (SMC). For a given set of 
proper motions, we investigate a range of orbital models 
based on different parametrizations. 

We assume that the Milky Way influences the orbits of 
the MCs through a fixed gravitational potential of three 
components: a central bulge, a disk, and an extended 
dark matter halo. For the halo, which is the dominant 
component of the MW at the large distances of the MCs, 
we have adopted an NFW density distribution (Navarro 
et al. 1996) suggested from CDM simulations: 



Po 



(cr/R vir )(l + cr/R vir y 



(1) 



where r is the spherical radius, po is the characteristic 
density, c is the concentration parameter, and i? v ; r is the 



virial radius. The total mass within r 
the virial mass and is given by 



i?vir is called 



Mvb - 4^ (iWc) 3 (ln(l + c) - c/(l + c)). (2) 

For halos with a mass comparable to that of the MW, 
typical values for the concentration parameter fall in the 
range c = 10 — 17 (Klypin et al. 2002), and we choose 
c — 12 in the present study. We consider Af v ; r and i? v ; r 
to be free parameters, allowing us to investigate different 
models for the orbital evolution of the MCs within the 



halo of the MW. In particular, we searched for halos 
within the range M vir = 1 - 2 X 10 12 Mq. 

The disk of the Milky Way is represented by a 
Miyamoto-Nagai (1975) potential 



GM A 



disk 



\Jr 2 + (a + Vz 2 +b 2 f 



(3) 



where is the total mass of the disk, a and b are scale 
parameters that control the radial and vertical extent of 
the disk, respectively, and R — \J x 2 + y 2 . 

We adopt a spherical Hernquist (1990) model for the 
potential of the Galactic bulge, 



^bulgc — — - 



GMh 



c b 



(4) 



where Mb and Cb are the total mass and the scale length 
of the bulge, respectively. 

In the present study, the above parameters for the disk 
and bulge are fixed to the following values (see also Ta- 
ble [l]): M d = 5.0x 10 10 M , a = 3.5 kpc, and b = 0.35 
kpc for the disk; M b = 0.5 x 10 10 M Q and c b = 0.7 kpc 
for the bulge (Binney & Tremaine 2008). That is, we did 
not change these values in our parameter space search for 
the best tidal model. Even though the disk and bulge do 
not influence the orbits of the MCs as strongly as the 
dark matter halo, the bulge and disk are nevertheless 
important in determining the circular velocity V4i r of the 
Milky Way (see Fig 2j. We have investigated models 
having the IAU standard value V c \ r = 220 km s _1 (Kerr 



Tidal Origin of the Magellanic Stream 



■5 



& Lynden-Bell 1986), and we have also considered a per- 
haps more realistic value of V^r = 240 km s -1 (Reid et 
al. 2009; Reid & Brunthaler 2004; Sirko et al. 2004). 
By varying the halo parameters M v j r and i? v ir and con- 
sidering the added affect of the bulge and disk, we can 
construct a variety of plausible MW models for our in- 
vestigation of MC orbits. 

The parameter values of our best model are summa- 
rized in Table [T] Our chosen NFW halo has M vir = 130 
x 1O 1O M and i? V ir = 175 kpc, which when added to 
the bulge and disk gives a total mass of 1.73 x 10 12 M Q 
within r = 300 kpc of the Milky Way. The circular ve- 
locity is Vcii — 240 km s . In the following sections 
we will describe the results of this model, but first we 
briefly mention an "alternate model" for the Milky Way 
in which the evolution of the MCs is remarkably similar. 
In particular, the results of sections [3] and [4] are largely 
unchanged when we substitute M vir =190 x 10 10 M and 
R vir =269 kpc for the NFW halo. In this model, the to- 
tal mass of the Milky Way within 300 kpc is 2.06 x 10 12 
M Q and the circular velocity at the solar radius is 220 
km s _1 . The predictions of this alternate model may 
be so similar in part because the total mass within the 
current location of the MCs r < 55 kpc is identical to 
that of the adopted model, 6.7 x 10 11 M . The rotation 
curves of the adopted and alternate models are given in 
Fig [2] The alternate model highlights the important 
fact that our chosen parameterization of the Milky Way 
is not unique and in fact is quite flexible when consider- 
ing the orbital evolution of the MCs. Further details of 
the alternate model are given in the Appendix. 

For the purposes of orbit integration, the LMC and 
SMC are each assumed to have a Plummer potential of 
the form 



$mc(0 



-GM n 



(5) 



where M mc is the total mass of either the LMC or SMC 
(MC), r is the distance from the center of mass, and a mc 
is the scale radius of the MC. In our parameter space 
search, we fixed the values of the scale radii, but we al- 
lowed the masses to vary as free parameters. Table [T] 
gives the adopted values. 

Dynamical friction is assumed to operate separately 
on the orbits of the LMC and SMC as they pass through 
the Milky Way's dark matter halo. We account for this 
effect by adopting the Chandrasekhar formula (Binncy 
& Tremaine 2008): 



■ffric 



47rG 2 M 2 c ln(A) /9dm (r) 



2A 



erf (A) = exp(-X 2 )]- 

u~ y7T 

(6) 

where M mc is the mass of either the LMC or the SMC 
(MC), v is the velocity of the MC, and A = v/(V2a), 
where a is the one-dimensional velocity dispersion of the 
adopted dark matter halo. As the MC moves through the 
halo, we calculate a at a given position using the analyt- 
ical approximation derived by Zentner & Bullock (2003). 
We adopt a reasonable value of 3.0 for the Coulomb log- 
arithm A (GN96). 

We carry out our simulations in a galactocentric frame 
such that the center of the Milky Way is always set to 
be (X,Y,Z) — (0,0,0). The initial position vectors and 
space velocities of the MCs are calculated from our choice 
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Fig. 2. — The rotation curve of the Milky Way, shown for its 
individual components (Bulge, Disk, NFW Halo) and the sum of 
the components (Total). A range of NFW halos were investigated, 
represented by the grey shaded regions. The adopted Milky Way 
model is given by the upper bound to the shaded regions (solid 
line; NFW halo M vir = 1.30 X 10 12 M , R vil =175 kpc), while an 
alternate model (see text) is indicated by the lower bound of the 
shaded regions (dashed line; NFW halo M vir = 1.90 X 10 12 M Q , 
-R v i r =269 kpc). The vertical dotted line at R = 8.5 kpc indicates 
the radius of the Sun, and its intersection with the Total rotation 
curve accordingly gives the circular velocity V c i r value. 

of free parameters (such as proper motions and V c ir] see 
Table fl] for adopted values) according to the method de- 
scribea in section 9 of van der Marel et al. (2002). We 
then integrate the appropriate equations of motion from 
the present epoch to 5 Gyr in the past using the method 
devised by MF80. By varying our free parameters, we 
are able to investigate a large number ~ 8 x 10 6 of ob- 
servationally constrained orbital models for the MCs. 

A number of studies suggest that the LMC and SMC 
became dynamically coupled only recently (Harris & 
Zaritsky 2009; Bekki et al. 2004), and in our previous 
work we explored the idea that the first strong inter- 
action during this coupling was responsible for creating 
the Magellanic Stream (DBlla). Guided by these re- 
sults, we restrict our attention to the MC orbital models 
that satisfy an assumed set of physical conditions for the 
formation of the Magellanic Stream. Namely, those as- 
sumed conditions are: (i) the mutual separation of the 
MCs has been less than 30 kpc for the majority of the 
last 2-3 Gyr, allowing them to dynamically interact, (ii) 
they strongly interact only twice, at t = —2.0 ±0.25 Gyr 
and t = -0.4 ± 0.25 Gyr (Harris & Zaritsky 2009), and 
(iii) their mutual separation cannot be too small (> 5 
kpc) at each epoch of strong interaction. These three 
requirements enable us to reduce the number of possible 
orbital models to ~ 10 5 . The models within this reduced 
set are further explored using test particle simulations in 
the next stage of investigation. 

2.2. Test Particle Simulations 

Using the set of orbits derived in the previous step, 
we run a large number of test particle simulations (e.g., 
see DB11) to determine which models are best able to re- 
produce the observed structure of the Magellanic Stream. 
Because our focus is mainly on the tidal evolution of the 
SMC, we represent the SMC by a disk of test particles 
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TABLE 2 
SMC N-body Parameters 





Halo 


Disk 


Spheroid 


Model 1 


AT dm =500,000 
M dm =1.36 x 10 9 M Q 
Rdm=7.5 kpc 


AT d =400,000 
M d =1.36 x 10 9 M Q 
R d =5.0 kpc 
6 d = -45° ; <f> d = 230° 


A r s P h=100,000 
M aph = 0.27 x 10 9 M 
R sp h=7.5 kpc 


Model 2* 


R dm =5.0 kpc 


i? d =5.0 kpc 


R sph =5.0 kpc 


Model 3* 


Rdm=5.0 kpc 


i? d =5.0 kpc 


R sph =2.5 kpc 



* Models 2 and 3 have the same parameters as Model 1 except for those noted. 
1 For brevity, variables are written without the subscript 'SMC as in the text. 
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Fig. 3. — Orbital separations between the LMC and MW (red), 
SMC and MW (green), and LMC and SMC (black). Solid lines 
indicate the time window in which the N-body models are evolved, 
t > —3.37 Gyr. Dotted lines indicate the continuation of the orbit 
integrations to 5 Gyr in the past. 

and leave the LMC as a fixed Plummer potential, as de- 
scribed previously. The SMC disk is assumed to have an 
exponential profile with a truncation radius of i?d,SMCj 
scale length of i?o,SMC = 0.2i?d,SMC, and a scale height 
of i?o,SMC = 0.02i? d sMC- Each test particle in the disk 
is given an initial circular velocity according to the rota- 
tion curve derived from the adopted Plummer potential 
for the SMC. In order to allow the investigation of a large 
number of models, we use a coarse resolution of 2000 to- 
tal particles in the test particle disk. This number is 
however more than enough to determine which orbital 
models provide promising tidal stripping scenarios for 
the MS. 

This phase of the study resembles the work of Ruzicka 
et al. (2010), who also perform a large suite of test par- 
ticle models in search of a formation scenario for the MS. 
However, the presently explored range of proper motions 
is much larger (Fig[T]) than that of Ruzicka et al. (2010), 
who restrict their attention to the la region of K06 and 
Piatek et al. (2008). Another key difference is that Ruz- 
icka et al. (2010) use an automated genetic algorithm 
to determine which models in their parameter space are 
best. In contrast, we use an ad hoc set of filters to intu- 
itively search through the candidate test particle models. 
In particular, we search for models that can satisfy the 
following requirements: (i) about 50% of the stripped 
SMC particles are located along the observed MS and 



(ii) more than 30% yet less than 70% of the SMC parti- 
cles are stripped. Inspecting the models by eye verified 
that these criteria could serve as good filters, but they 
are by no means unique or optimal. Our intent is to 
simply show how we arrived at our preferred model. 

Imposing the above criteria enabled us to further nar- 
row the parameter space, as only ~10% of the test parti- 
cle models pass these conditions. Of these models, most 
fall into one of several "families" which exhibit similar 
tidal stripping of the SMC disk. We take the best can- 
didates from these families and pass them to higher res- 
olution N-body simulations to provide more faithful dy- 
namical models. 

2.3. N-body Simulations 

In this final phase of the investigation, we represent 
the SMC by a system of N-body particles whereas the 
same fixed Plummer potential is used for the LMC. We 
search for the best few models which can explain the 
observed properties of the MS, LA, Bridge, and SMC at 
the present day. 

Almost all previous investigations of the formation of 
the MS adopted a bulgeless disk ("pure disk") for the 
SMC (e.g., GN96, YN03, C06, B10, and DBlla,b), and 
we accordingly adopt a pure disk model as a first step. 
However, recent observations have suggested that the 
older stellar populations of the SMC are supported by 
their velocity dispersion rather than rotation, implying 
that the SMC may possess a spheroid population of old 
stars (Harris & Zartisky 2006). Considering this observa- 
tional result, our final investigation involves "disk-plus- 
spheroid" models in which the baryonic component of 
the SMC is composed of a disk and a central spheroid. 
We follow a two-step process: first, we run medium res- 
olution N-body simulations representing the SMC by a 
self-gravitating disk, and we thereby identify the best 
tidal models from the previous set of test particle simu- 
lations; and second, we represent the SMC by a multi- 
component disk-plus-spheroid system and run a final set 
of high resolution simulations for our adopted model. In 
the present work, we ignore drag forces induced by the 
hot halo of the Milky Way (e.g., Diaz & Bekki 2011b) and 
focus instead on the essential gravitational dynamics of 
the system. 

Numerical computations were carried out at the Uni- 
versity of Western Australia on (i) the latest version 
of GRAPE (GRavity PipE, GRAPE-DR), which is 
a special-purpose computer for gravitational dynamics 
(Sugimoto et al. 1990), and (ii) a Core i5 desktop com- 
puter system with a GPU card (NVIDIA GTX580) im- 
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plementing the CUDA G5/G6 software package for cal- 
culations of gravitational dynamics. 

The time integration of the equations of motion is per- 
formed by using a second-order leapfrog method with a 
time step interval of ~ 0.02i dyn , where i dyn is the dy- 
namical time scale of the SMC. To save on computa- 
tional cost, we do not simulate over the full 5 Gyr of 
the computed orbit. Instead, we begin our simulations 
1 Gyr or so prior to the first strong interaction between 
the SMC and LMC. For our adopted model, the sim- 
ulated time window was 3.37 Gyr (see Appendix for 
discussion). The total number of particles for the first 
set of models (medium resolution, pure disk models) is 
X=200,000, split evenly between the disk and halo com- 
ponents. The total number of particles in the second 
set of models (high resolution, disk-plus-spheroid mod- 
els) is N = 10 6 , where A r dm =500,000 for the dark matter 
halo, N d =400,000 for the disk, and iV sph =100,000 for 
the spheroid. The gravitational softening lengths e for 
all components are set to be 72 pc. 

2.3.1. Pure Disk Models 

In this intermediate step we represent the SMC as a 
bulgeless disk galaxy embedded in a massive dark matter 
halo. The total mass and the size of the dark matter halo 
(disk) is M dmjSMC (M diS Mc) and i? dm ,s M c (#d,SMc), re- 
spectively We adopt an NFW profile for the density dis- 
tribution of the halo where the scale length r Sj sMC is set 
to be 0.5i? dj sivic- The mass ratio of dark matter to total 
SMC mass is fixed at / dm ,SMC = 0.5, which is consistent 
with the observed HI rotation curve for the SMC (e.g., 
Bekki & Stanimirovic 2009). Although we did not ex- 
plore other values for / d m,SMC, we should point out that 
a wide range of values could be consistent with the HI 
rotation curve depending on stellar mass-to-light ratios. 

The radial (R) density profile of the SMC disk is as- 
sumed to be proportional to exp(— R/ Rq,smc) with scale 
length Rq,smc — 0.2i? dj sMC> an d the vertical (Z) density 
profile is assumed to be proportional to sech 2 (Z/Zo.sMc) 
with scale length Zq^smc = 0.04i? d .sMC- I n addition to 
the rotational velocity induced by the gravitational field 
of the disk and halo, the disk is assigned both radial and 
azimuthal velocity dispersions according to the epicyclic 
theory with Toomre's parameter Q = 1.5. The vertical 
velocity dispersion at a given radius is set to be 0.5 times 
as large as the radial velocity dispersion at that point. 

Following the methodology of GN96, the spin of the 
SMC disk is specified by two angles: d , which is the an- 
gle between the Z-axis of the Milky Way and the angular 
momentum vector of the SMC disk, and <f>&, which is the 
azimuthal angle. In particular, </> d is measured from the 
X-axis to the projection of the angular momentum vec- 
tor onto the X-Y plane of the Milky Way (see Fig 1 of 
GN96 for schematic). We investigated the full range of 
values for # d and ^ d for each of our best models, and 
we discuss the range of preferred values in the Appendix. 
After analyzing the full set of pure disk models, we ended 
up with a final set of two best models, our adopted model 
and an alternate model. Although these two models are 
remarkably similar in terms of how well they reproduce 
the MS, LA, and Bridge, we will focus mostly on the 
adopted model and leave the alternate model to the Ap- 
pendix. Table [I] gives the final set of parameters, and 



Fig [3] gives the orbital separations for the adopted model 
over 5 Gyr. 

2.3.2. Disk-plus-Spheroid Models 

Whereas the numerical investigations up to this point 
were solely for the purpose of finding the best model, the 
N-body simulations described here are presented as our 
main results in sections [3] and [4] We represent the SMC 
as a mult i- component system composed of an NFW dark 
matter halo, an exponential disk, and a central spheroid, 
similar to the models of Bekki & Chiba (2009). The 
N-body prescriptions for the dark matter halo and the 
exponential disk are the same as those described in the 
pure disk model. The spheroid is represented by a Plum- 
mer model with mass M sp h,SMCj scale length a sp h,SMCj 
and truncation radius -R S ph,SMC- We make the simple as- 
sumption that the spheroid pertains to an old population 
of "stars" and that the disk (particularly its outer parts) 
pertains to "gas." We therefore investigate whether the 
particles stripped from the disk are able to reproduce the 
observed properties of the MS. 

Although Nidever et al. (2011) have recently inves- 
tigated the outer structure of metal-poor giant branch 
stars in the SMC, it remains observationally unclear how 
stellar populations of different ages and metallicities dis- 
tribute within the tidal radius of the SMC. We therefore 
assume that i? S ph,SMC an d a S ph,SMC are free parameters. 
Likewise, we consider the mass ratio of the spheroid to 
the disk / sp h,SMC to be a freely adjustable parameter. In 
our investigation of the formation processes of the MS, 
we find that reasonable variations in the mass ratio and 
scale length of the spheroid do not strongly affect its tidal 
evolution, and we accordingly adopt / sp h,SMC = 0.2 and 
& S ph,SMC = 0.2i? sp h,SMC- However, changes in the trun- 
cation size of the spheroid strongly alter the tidal evolu- 
tion, as described in section [4] Accordingly, we present 
the results of three high resolution models where the size 
of the spheroid is varied: an extended model ("Model 1" 
-Rsph.SMC = 7.5 kpc), an intermediate model ("Model 2" 
-Rsph.SMC = 5 kpc), and a compact model ("Model 3" 
-Rsph.SMC = 2.5 kpc). Parameter choices for our best 
disk-plus-spheroid models are summarized in Table [2] 

Despite the simplicity of the present study (e.g., ne- 
glect of gas dynamics, star formation, chemical evolu- 
tion), this work is nevertheless an important step toward 
constructing a self-consistent model for the formation of 
the MS and the evolution of the Magellanic system as a 
whole. First, we have run a large number of models and 
found a promising set of parameters (including a set of 
proper motions for the MCs) which may well inform fu- 
ture modeling. And second, as detailed in the following 
sections, we are able to explain many of the properties of 
the Magellanic system and MS using collisionless dynam- 
ics alone. Where the model is unsuccessful, we consider 
the inclusion of gas dynamics to be a salie nt o ption for 



future improvements (see Discussion section 5.2). In par- 
ticular, insights will surely be provided by modeling the 
hydrodynamical interaction of MS gas clumps with the 
hot halo of the Milky Way (e.g., Bland-Hawthorn et al. 
2007). 

3. RESULTS I - THE SMC DISK 

Here we describe our results concerning the SMC disk 
of Model 1. See Tables [l] and [2] for parameter sets. 
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Fig. 4. — The tidal evolution of the SMC disk. Each panel shows the distribution of SMC disk particles in the y-z plane at a specific 
time. Also displayed are the orbit trails of the SMC (green) and LMC (red). The Milky Way is shown schematically as an edge-on disk 
with radius 20 kpc, and the position of the Sun at (-8.5, 0, 0) kpc is indicated by the blue circle. Because of the chosen projection in the 
y-z plane, the location of the Sun appears to coincide with the Galactic center (0, 0, 0). In each of the first two panels an inset is provided 
which shows the face-on state of the SMC disk. At t = —2.48 Gyr (top left panel and inset) the SMC disk has attained spiral structure 
due to weak tidal interaction with the LMC and Milky Way. At t = —1.92 Gyr (top middle panel and inset) the SMC disk is disrupted by 
a close passage of the LMC, causing two tidal arms (A and B) to be stripped away. At t = —1.30 Gyr (top right panel) arm A lags behind 
the SMC orbit while arm B is engulfed by the LMC. Eventually, arm B forms a ring structure in the LMC (t = —0.68 Gyr, bottom left 
panel). The leading and trailing debris elongate as the Milky Way pericenter is approached (t = —0.34 Gyr, bottom middle panel), and by 
t = Gyr (bottom right panel) the Magellanic Stream (MS) and Leading Arm (LA) are fully formed. Not pictured is the formation of the 
Magellanic Bridge and Counter-Bridge, which happens at t = —0.26 Gyr as described in the text. 



3.1. Tidal Evolution of the SMC Disk 

In our adopted model, the LMC and SMC have suffered 
only two close passages of one another, at t = —1.97 Gyr 
(the first interaction) and t = —0.26 Gyr (the second in- 
teraction) as seen in Fig[3j Previous to these encounters, 
the LMC and SMC were separated by large distances (up 
to ~ 120 kpc at t = — 5 Gyr) as they orbited indepen- 
dently within the virial radius of the Milky Way halo. In 
this sense, the formation of a strong binary state is a re- 
cent phenomenon for the MCs, and this scenario provides 
a framework for understanding much of the dynamics. 
For example, the strong tidal forces suffered by the SMC 
during the first interaction are responsible for stripping 
away much of its disk. Figure H] shows the time evolution 
of the SMC disk along with the orbital evolution of the 
MCs about the MW, and a snapshot of the disruption of 
the SMC is given in the top middle panel. 

Previous to this epoch of strong interaction, the SMC 
is subject to a combination of weak tidal fields from the 
LMC and MW (e.g., at t < -2 Gyr, Fig [j§. These 
weak tidal forces are sufficient to induce a bar and spiral 
arms within the SMC (Fig |4j top left panel) , but they 
are not strong enough to disrupt the disk. In particular, 
even though the MW pericenter at t = —2.3 Gyr pre- 
cipitates the formation of the strong LMC-SMC binary 



state, the MW pericenter itself does not actively par- 
ticipate in the removal of material from the SMC. The 
disk remains intact up until the first close passage of the 
LMC at t = —1.97 Gyr, when the mutual separation be- 
tween the MCs reaches a minimum of 6.0 kpc. At this 
point, two tidal arms ("A" and "B", top middle panel of 
FigH begin to separate from the main body of the SMC. 
These two tidal arms follow quite different evolutionary 
paths: arm A falls into orbit behind the SMC and slowly 
elongates to form the Magellanic Stream, whereas arm 
B stretches ahead of the orbit and is largely engulfed 
by the LMC (Fig|4j top right an d bott om left panels). 
As explored in the L)iscussion Sec |5.3.l| this transfer of 
SMC disk material to the LMC would fuel star forma- 
tion via the infall of metal-poor gas. At t = —0.68 Gyr 
(Fig [I} bottom left panel) the captured mass within the 
LMCresembles a polar ring, suggestive of the distribu- 
tion of carbon stars in the present-day LMC (Kunkel et 
al. 1997a). 

By t = —0.34 Gyr, the MW has accelerated a fraction 
of the stripped material into the leading direction of the 
orbit (Fig HI bottom middle panel). Though not imme- 
diately clear from Fig [4] the leading material is collected 
into two separate branches, one that originates in tidal 
arm A and the other in B. In Sec 13.41 we will indicate 
that these branches are spatially and kinematically dis- 
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Fig. 5. — (Left panel) The final distribution of disk particles projected onto the sky. Column densities are computed as a smoothed 
intensity map of the simulated mass distribution. The coordinates are galactic longitude I (straight lines) and latitude b (concentric circles) 
shown in a ZEA projection centered on the South Galactic Pole. The orbit trails of the SMC (green) and LMC (red) are also shown, 
as are their current locations (circles). (Right panel) The observed HI column densities of the Magellanic system, combined from data 
presented in Putman et al. (2003a) and Nidever et al. (2010). The densest regions (yellow) are the SMC (I = 302.8°, b = -44.6°) and 
LMC (I = 280.5°, b = —32.5°). The Magellanic Stream is the prominent trail of gas that splits into two filaments parallel to / = 270° and 
I = 90°. Other notable features include the Magellanic Bridge, which extends between the LMC and SMC, and the many branches of the 
Leading Arm at the bottom of the panel. Column densities are represented on a logarithmic scale as shown at the base of the figure. 

tinct. At t = -0.26 Gyr, the LMC and SMC undergo 
their second strong interaction, reaching a minimum sep- 
aration of 6.6 kpc. This second close passage strips an 
additional two tidal arms from the SMC disk, but only 
one of these structures corresponds to an observable HI 
feature of the Magellanic system, namely the Magellanic 
Bridge. The other structure, which we call the "Counter- 
Bridge," has perhaps eluded identification by extending 
almost directly behind the SMC along the line of sight. 
The three-dimensional natu re of the Bridge and Counter- 
Bridge is considered in Sec |3.5| At t = Gyr, the MCs 
have arrived to their present-day positions and velocities, 
and the simulated particle distribution has reached its fi- 
nal state. The locations of the simulated Leading Arm 
and Magellanic Stream are noted in the bottom right 
panel of Fig [I] 



3.2. Projected On-sky Distribution 

To make a meaningful comparison with observations, 
we must project the final distribution of particles against 
the sky and create a map of their density. The left hand 
panel of Fig [5] gives the simulated column density distri- 
bution in a ZEA projection centered on the south galactic 
pole. The coordinates are galactic longitude (I) and lati- 
tude (6). The orbit trails of the SMC and LMC are pro- 
vided (green and red, respectively), as are their current 
locations (circles). The right hand panel gives the ob- 
served HI column density map of the Magellanic system, 
where column densities are plotted on the same logarith- 
mic scale as the simulation data. Although the observed 
LMC is rich in HI (right panel), our simulated column 
density map does not display a high density region at the 
location of the LMC. This is simply due to our model as- 
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Fig. 6. — The evolution of the Magellanic Stream, shown as a 
time series starting from the inset (top left, t = —1.92 Gyr) and 
proceeding to the left panel (f = —0.90 Gyr), then the middle 
(t = —0.49 Gyr), and finally the right (t = Gyr). The coloring is 
chosen at t = —1.92 Gyr (inset), where the particles in tidal arm 
A (see Figure El are colored according to their distance along the 
tidal arm, andall other particles are black. This coloring is pre- 
served at each time step as the tidal arm separates from the SMC 
(left and middle panels) and finally coincides with the present day 
Magellanic Stream (right panel). See text for an analysis including 
a description of the dynamics. The projection used in the right 
panel (t = Gyr) is an on-sky projection in galactic coordinates 
(Z, b), and for ease of reference the same projection is used for the 
middle and left panels. The SMC and LMC are represented as cir- 
cles (green and red, respectively). The inset is taken from Figure[4] 
albeit without color. 

sumptions, because the LMC is represented by a fixed 
potential rather than a live ensemble of particles. The 
spatial extent of HI in the SMC and Bridge appears to 
be larger in the simulation than in observation, which is 
also the case in the collisionless simulations of Connors 
et al. (2006) . This discrepancy could be mitigated by in- 
cluding gaseous pressure from the Milky Way's hot halo, 
which would confine the gas. 

Judging from Fig [5l the simulation provides a good 
reconstruction of the Bridge and MS in terms of on-sky 
location and morphological structure. The leading ma- 
terial, however, does not provide a good matc h to the 
observed location of the LA. We argue in Sec |5.2| that 
the inclusion of ram pressure from the Milky Way s hot 
halo may be able to resolve this discrepancy. The tip 
of the MS also poses a problem, because the simulation 
over-predicts the density as compared to observations, 
particularly for b > —40°. It would appear from Fig [5] 
that the observed MS does not even extend into this dis- 
tant region, but Nidever et al. (2010) have shown that 
MS gas does indeed exist here, albeit at low column den- 



sities. As compared to the observations of Nidever et al. 
(2010), our simulation correctly predicts the total length 
of the stripped Magellanic gas (~ 200° from MS tip to 
the end of the LA), and it also correct ly p redicts the 
kinematical structure of the MS tip (Sec |3.4[ ). 

The triangle-shaped structure at (I, b) = (315°, —50°) 
(Fig [5j right panel) is reproduced reasonably well in 
our model, despite having a slightly different orienta- 
tion. The structure has not yet been given a unique 
identity in either the observational or theoretical litera- 
ture: no simulation has yet addressed its formation, and 
observers have associated it with the gas clouds at the 
base of the MS (e.g., Putman et al. 2003a; Bruns et al. 
2005, who use the collective label "Interface Region"). 
We find that the origin of this structure is quite different 
to that of the MS, and we therefore propose it be named 
to reflect its unique identity. In the present work we will 
use the name "Magellanic Horn," owing to its tapered 
appearance. Whereas the MS is a product of the elon- 
gation of tidal arm A, the Horn is a remnant of tidal 
arm B. When portions of tidal arm B are accreted onto 
the LMC, much of the additional material is left in orbit 
around the LMC periphery. Some of this peripheral ma- 
terial is slung away as the LMC and SMC plunge toward 
one another at t — —0.26 Gyr, becoming unbound from 
the MCs. The material that is left behind in the trail- 
ing direction of the orbit gradually elongates to form the 
Horn. 

Perhaps the most compelling aspect of the simulation 
is the reproduction of MS bifurcation. As seen in the 
right panel of Fig[5j the base of the MS bifurcates into 
two parallel and distinct filaments, and the location of 
these filaments is well reproduced in the simulation. An 
important observational feature of the filaments is that 
they appear to cross at several locations along the MS 
(Putman et al. 2003a). The location of the first cross- 
ing point at (I, b) « (45°, —80°) is clearly exhibited in 
the simulated MS, providing a convincing morphological 
match. Considering the strong reproduction of the MS 
filaments, it will be instructive to describe their origin 
in detail. We do so in the following section, and we 
co mpar e with the findings of Nidever et al. (2008) in 
Sec IQl 

3.3. The Origin of MS Bifurcation 

Figure [6] color-codes the time-evolution of the MS. To 
construct the color-code, we first fit a logarithmic spiral 
to tidal arm A at t = —1.92 Gyr (Figpj inset of top 
middle panel), and then we assign color to the particles 
based on their projected distance along the arm (Fig [6j 
inset). The colors run across a portion of the HSV spec- 
trum beginning at the base of tidal arm A (purple and 
blue) continuing to its middle (cyan and green) and end- 
ing at its tip (yellow and orange). The particles retain 
their color throughout the evolution, allowing us to trace 
the origin of structures as the MS is formed. 

As seen in the right panel (t — Gyr) of Fig [6j the 
two filaments of the MS originate from opposite ends of 
tidal arm A. The right filament originates from the base 
of arm A (purple and dark blue) whereas the left fila- 
ment comes from its middle and end (green and yellow, 
respectively). In addition, the colors that fall in between 
these extremes are found in the MS tip (blue, cyan). Ac- 
cordingly, we can think of the MS as being a wrapping 
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of tidal arm A: (i) the base of arm A becomes the base 
of the right MS filament, (ii) arm A is traced upward 
through the right filament toward the MS tip, and (iii) 
the middle of arm A is traced downward through the left 
MS filament. Though it appears faint in Fig|6j the tip 
of arm A (orange) continues to trace downward through 
the left filament into the leading region. The st ructu re 
of the leading material is made more clear in Sec |3.4| 

More than simply tracing up and back down the MS, 
the material of tidal arm A also crosses within the body 
of the MS. Consider the left panel of Fig [6] (f = -0.90 
Gyr) , in which the material from the base of arm A (pur- 
ple and blue) is to the left of that from the end of arm 
A (green and yellow). In the middle panel (t = —0.49 
Gyr) , the material from the base overlaps with that from 
the end of arm A, and finally in the right panel (t = 
Gyr) , the material from the base is to the right of that 
from the end. In other words, the relative on-sky orien- 
tation of the MS filaments has reversed within the last 
~1 Gyr. This scenario provides a physical interpretation 
for the crossing point at (I, b) « (45°, —80°): it is a re- 
gion of physical overlap, similar to the crossover point of 
a ribbon that has been folded across itself. 

But what dynamical processes are responsible? One 
important point of emphasis is that the tip of the MS 
(blue) does not come from the tip of tidal arm A (yellow- 
orange). This idea runs counter to the expectation from 
two-body tidal stripping models (e.g., Toomre & Toomre 
1972; Besla et al. 2010) in which the tips of tidal tails 
elongate away from their parent galaxies. The complex 
evolution of tidal arm A underlines the three-body na- 
ture of the formation of the MS. Tidal arm A is removed 
from the SMC disk due to a close passage of the LMC, 
but once it becomes unbound, its subsequent evolution 
is dominated by the MW. In particular, the approach to 
MW pericenter over the past ~1 Gyr is responsible for 
elongating tidal arm A and creating an identifiable MS 
tip, as seen in the time progression of Fig[6] Because they 
are created by distinct dynamical processes, there is no 
reason to expect that the tip of the MS should coincide 
with the tip of tidal arm A. 

The MW also accelerates much of the unbound mate- 
rial of tidal arm A into the leading direction of the orbit. 
In the left panel of FigJTj] (t = -0.90 Gyr), the entirety 
of arm A is on the trailing side of the orbit, but as MW 
pericenter is approached (middle and right panels), the 
tip of arm A (yellow-orange) is accelerated to the leading 
side. Notice that on its path to the leading side (as seen 
by comparing the three main panels of Fig |6| , the tip of 
arm A (yellow-orange) must cross from the right side to 
the left side of the purple/blue filament. The MW coaxes 
the adjoining material of arm A (i.e., the green/yellow 
filament) to follow a similar evolutionary path, causing it 
also to cross over. The material from the base of arm A 
(i.e., the purple/blue filament) is perhaps less susceptible 
to this effect because it originates closer to the SMC at 
the base of tidal arm A, and it is therefore more strongly 
anchored in place. 

3.4. Distance and Kinematics 

Fig [7] gives the distances and vlsr radial velocities of 
the simulation particles, and the color scheme of Fig [6] 
has been retained in order to permit the identification of 
various structures. The MS (colored) and Horn (black) 



are found on the trailingside of the orbit (Lms < —20°; 
right hand side of Fig [7]), and the leading material is 
found on the leading side of the orbit (Lms > 10°; left 
hand side of Fig [7]). It is clear that the leading mate- 
rial has two distinct branches, one that is sourced from 
tidal arm A (i.e., fro m its tip, colored yellow-orange as 
described in Sec 3.3), and the other from tidal arm B 
(colored black). A comparison of the top and bottom 
panels of Fig [7] shows that these two branches are both 
kinematically and spatially distinct, despite being over- 
lapped in the on-sky map of Fig[5j The two branches are 
not separated on the sky simply because of a projection 
effect. 

The observed HI kinematics of the Magellanic system 
are given as grey contours in the background of Fig [7J 
lower panel (Nidever et al. 2010). The correspondence 
between simulation and observation is good for the MS. 
The main body of the MS correctly predicts the observed 
radial velocity profile, and the two MS filaments pos- 
sess an essentially identical (i.e., overlapping) kinemati- 
cal structure. The tip of the MS shows a distinct velocity 
inflection that is echoed in the observations of Nidever 
et al. (2010). This fact provides a good justification for 
the material at the tip of the simulated MS: even though 
it has the wrong density (see Fig [5]), it has the correct 
extent and the correct kinematics. The radial velocity 
profile of the Horn is slightly offset to that of the MS, 
but this is to be expected because the Horn originates 
in tidal arm B along a very different evolutionary track 
(see Sec 3.2 ) as compared to the MS, which originates in 
arm A. 

The correspondence between simulation and observa- 
tion is comparatively poor for the LA. The shape of the 
velcocity profile is generally correct for both branches 
of leading material, but it is offset by ~ 150 km s" 1 
from observations. The heliocentric distances (FigJTJ up- 
per panel) also suggest a discrepancy for the LA region. 
Whereas the simulation predicts that the LA should lie 
at distances of ~ 50 kpc to ~ 100 kpc, McClure-Griffiths 
et al. (2008) have provided a much smaller distance of 
~20 kpc. The uncertainty of this data point is consid- 
erable, because it is a kinematic distance derived from a 
supposed interaction between the LA and the HI disk of 
the MW. Nevertheless, the large discrepancy in distances 
for the LA, coupled with the failure to correctly predict 
the kinematics and on-sky location, suggests that we are 
missing an es sent ial ingredient in the dynamics. As dis- 
cussed in Sec |5.2[ we consider this ingredient to be ram 
pressure from the MW hot halo. 

Observational data does not appear in the upper 
panel of Fig [7J because distances to HVCs cannot be di- 
rectly measured, barring the one data point of McClure- 
Griffiths et al. (2008) for the LA. At best, only upper and 
lower bounds can be obtained by carefully considering 
foreground/background absorption spectra (e.g., Wakker 
2001), but this has not been done extensively for the MS 
and LA, whose distance profiles are completely unknown 
(however, see Jin & Lynden-Bell 2008). Accordingly, it 
cannot yet be verified if the MS has a large radial extent 
as predicted in the simulation (Fig[7J upper panel). This 
large radial spread is quite curious, especially when we 
consider that the on-sky projection of the MS is filamen- 
tary. It may be more appropriate to describe the MS as 
radially extended sheets rather than confined filaments. 
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Fig. 7. — The heliocentric distances (top panel) and dlsr velocities (bottom panel) of the disk particles. The horizontal axis is Magellanic 
longitude Ljyis, which runs parallel to the Magellanic Stream and is defined by Nidever et al. (2008). The coloring of the particles is 
maintained from Figure m\ In each panel, the Leading Arm is to the left of the SMC (green circle) and LMC (red circle), whereas the 
Magellanic Stream is to fne right. In the bottom panel, observational data from Nidever et al. (2010) is represented by grey filled contours 
(at 1, 10" 10 , and 10~ 25 Kelvin deg). 



The responsible dynamical process is simple: at the same 
time that the MW elongates the MS on the sky over the 
past ~1 Gyr, radial tidal forces from the MW cause elon- 
gation along the line of sight. These processes cannot be 
decoupled, but perhaps the elongation is mitigated by 
hydrodynamical interactions such as cloud confinement 
and ablation (e.g., Bland-Hawthorn et al. 2007). 

3.5. Magellanic Bridge and Counter-Bridge 

It is clear that the first interaction between the MCs 
(at t = —1.97 Gyr) removes two tidal structures from 
the SMC disk (e.g., Fig[2|), but it is less clear that an 
additional two structures are removed during the sec- 
ond strong interaction (at t = —0.26 Gyr). At least one 
tidal structure obviously forms, because it fills the region 
in between the SMC and LMC in good agreement with 
the observed Magellanic Bridge. The other structure is 
less obvious, however, and here we call it the "Counter- 
Bridge" . Panel (a) of Figure [8] shows a close-up view of 
the simulated column density map of Fig [5] centered on 
the Magellanic Bridge. In this panel, the existence of 
the Counter-Bridge is traced only by a small region of 
high density around —Lms ~ 16°, —Bms ~ 6°. Pan- 
els (b) and (c) reveal the three-dimensional structure of 
the region by plotting the heliocentric distances of the 
simulation particles along each spatial direction. These 
panels reveal a dense and clearly defined tidal structure, 
i.e., the Counter-Brige, extending away from the SMC 



toward distances of up to d ~ 85 kpc. 

Despite the ~ 20 kpc extent of the Counter-Bridge, 
it remains propitiously aligned with the SMC along the 
line of sight and is therefore hidden in the on-sky col- 
umn density distribution. Perhaps for this reason alone 
the Counter-Bridge has eluded observational identifica- 
tion. However, a "loop" of HI on the northeast edge 
of the SMC may betray the existence of the Counter- 
Bridge, as suggested by Muller & Bekki (2007). Another 
possible observational analogue of the Counter-Bridge is 
the large line-of-sight depth of the SMC, which has been 
measured to be ~ 6 to ~ 12 kpc (Welch et al. 1987; 
Crowl et al. 2001; Lah et al. 2005), - 12 to - 16 kpc 
(Hatzidimitriou & Hawkins 1989; Gardiner & Hawkins 
1991; Subramanian & Subramaniam 2012), or up to ~ 20 
kpc (Mathewson et al. 1988; Groenewegen 2000). De- 
spite the disagreement in these observational measure- 
ments (see HZ06 for a clarifying discussion), it is clear 
than the SMC is significantly extended along the line of 
sight and particularly in the north-east quadrant. The 
present model provides a natural mechanism to explain 
this structure, namely that a projection effect along the 
line of sight causes the bound stellar populations of the 
SMC to be confused with the unbound stars tracing the 
tidally extended Counter-Bridge. 

The radial velocity signature of the Counter-Bridge is 
given in Fig [91 which plots ulsr against heliocentric dis- 
tance D for all particles that roughly overlap the SMC 
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Fig. 8. — The Magellanic Bridge and Counter-Bridge, (a) A 
close-up of the column density map of Figurepjlusing the same color 
scheme. The panel is centered on the Magellanic Bridge, which 
extends between the SMC (green circle) and LMC (red circle). 
The on-sky coordinates are Magellanic longitude Lms an d latitude 
Bms as defined by Nidever et al. (2008). Contours are drawn at 
the 0.96, 0.98, and 0.995 quantiles, outlining the Bridge and its 
tidal complement the Counter-Bridge, which is located at (Lms, 
Bms) ~ ( — 16°, — 6°). The heliocentric distances of the simulation 
particles are plotted against (b) Magellanic latitude Bms, an d ( c ) 
Magellanic longitude Lms- The coordinates of the SMC and LMC 
are shown as green and red circles, respectively. The presence of 
the Counter-Bridge becomes more obvious in panel (b) where it is 
located at top right, between Bms ~ — 4° and —16°, and in panel 
(c) where it is located at right, between Lms ~ — 14° and —20°. 
As seen in both panels, the Counter-Bridge extends up to 20 kpc 
away from the SMC largely along the line of sight. 

in the sky (-22° < L MS < -12°, -18° < B M s < -2°). 
There are three clear components, each residing at dis- 
tinct distances: particles bound to the SMC (59 kpc 
< D < 65 kpc), those in the Bridge (D < 59 kpc), 
and those in the Counter-Bridge (D > 65 kpc). Be- 
cause much of the Bridge is found elsewhere on the sky, 
the velocity profile for the Bridge is represented incom- 
pletely in Fig|9j On the other hand, the velocity profile of 
the Counter-Bridge is fully characterized in Fig[9j which 
shows a linear rise in velocity with distance. Despite this 
unique signature, the Counter-Bridge remains confused 
with the SMC in projection, because each exhibits the 
same range of radial velocities in roughly the same loca- 
tions on the sky. 

Accordingly, it may be difficult to observationally ver- 
ify/refute the existence of the Counter-Bridge without 
obtaining distance measurements. For this reason, HI 
observations are not as suitable to the task as stellar 
surveys. The presence of young stars in the Bridge 
(e.g., Irwin et al. 1985) may indicate that young stars 
trace the Counter-Bridge as well, owing to the fact that 
the Bridge and Counter-Bridge are complementary tidal 
structures with similar origins. The use of older stars 
as tracers of the Counter-Bridge may prove difficult be- 
cause intermediate-age stars have yet to be observed in 
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Fig. 9. — Velocity dlsr plotted as a function of heliocentric dis- 
tance D for particles originating within the disk of the SMC. The 
only particles plotted here are roughly coincident with the SMC 
as projected onto the sky (-22° < Lms < -12°, -18° < B M s < 
-2°). The location of the SMC (59 kpc < D < 65 kpc) and 
Counter-Bridge (D > 65 kpc; "CB") are labeled. 
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Fig. 10.— The v LS n velocity field of (a) the SMC disk of model 
1, and (b-d) the SMC spheroid of models 1-3, respectively. In each 
panel, the on-sky coordinates are Magellanic longitude Lms an d 
latitude Bms as defined by Nidever et al. (2008). The dashed 
circle has a radius of 3 kpc and is centered on the SMC center 
of mass. The velocity range has a lower bound of 100 km s 
(blue), an upper bound of 180 km s — 1 (red), and is centered on 
the ulsr velocity of the SMC (white). Only the disk (a) exhibits 
a large rotational amplitude, because the velocity gradient for the 
spheroids (b-d) is due to tidal effects (see text). 




Fig. 11. — The tidal evolution of the SMC spheroid. Each panel shows the distribution of SMC spheroid particles of models 1 (cyan), 2 
(black), and 3 (magenta) in the y-z plane at a specific time. Also displayed are the orbit trails of the SMC (green) and LMC (red). The 
Milky Way is shown schematically as an edge-on disk with radius 20 kpc, and the position of the Sun at (-8.5, 0, 0) kpc is indicated by 
the blue circle. Because of the chosen projection in the y-z plane, the location of the Sun appears to coincide with the Galactic center 
(0, 0, 0). In each of the first two panels an inset is provided which shows the face-on state of the SMC spheroid. At t = —2.48 Gyr (top left 
panel and inset) the SMC spheroid is only weakly distressed by tidal forces, and the key differences between model 1 (extended), model 2 
(intermediate), and model 3 (compact) are clear. The close passage of the LMC at t = —1.92 Gyr (top middle panel and inset) causes a 
number of particles to be stripped away, particularly for models 1 and 2. A portion of these stripped particles are engulfed by the LMC 
(t = —1.30 Gyr, top right), while others orbit freely in leading and trailing streams (t = —0.68 Gyr, bottom left). A second close passage at 
t = —0.34 Gyr (bottom middle) pulls more spheroid particles from the SMC. At the present day (t = Gyr, bottom right) stellar analogues 
for the Magellanic Bridge, Leading Arm, and Magellanic Stream have been formed. 



the Bridge (Harris 2007) . Regardless of the tracer pop- 
ulation used, establishing a distance-velocity correlation 
akin to Fig [9] would provide a strong test of the existence 
of the Counter-Bridge. 

3.6. SMC velocity profile 



Panel (a) of Fig 10 gives the vlsr velocity field of the 
SMC disk at the present day. There is a strong sign 
of rotation, consistent with the observed HI kinematics 
of the SMC (Stanimirovic et al. 2004). Rotation is 
not consistent, however, with the large population 
of red giant stars observed within the SMC (Harris & 
Zaritsky 2006). This old stellar population appears to be 
supported instead by its velocity dispersion, suggestive 
of a spheroidal morphology. We are therefore motivated 
to investigate the tidal evolution of such a spheroid in 
the following section. 



4. RESULTS II - THE SMC SPHEROID 

Here we describe our results concerning the SMC 
spheroid of Models 1, 2, and 3 (see Tables |l] and [2] for 
parameter sets). Because little is known about the evo- 
lution of the spheroid, we present multiple parameteri- 
zations as a comparative study. The essential difference 



between these spheroids is their initial size (i?sph.SMc): 
Model 1 is extended (=7.5 kpc), Model 2 is intermediate 
(=5.0 kpc), and Model 3 is compact (=2.5 kpc). 

4.1. Tidal evolution of the SMC spheroid 

The SMC sp hero id evolves under the same tidal forces 
as the disk (Sec 3.1 1, and its evolution is of course coupled 
to that of the disk through the N-body dynamics. Fig[TT1 
shows the tidal evolution of each spheroid model along 
with the orbits of the MCs. The particles belonging to 
models 1, 2, and 3 are colored cyan, black, and magenta, 
respectively. We showed previously that the SMC disk 
was disrupted following the first strong interaction be- 
tween the MCs at t — —1.97 Gyr, and we see from Fig [TT| 
(top middle panel) that the same forces are responsible 
for disrupting the spheroid. Only models 1 and 2 suffer 
significant mass loss, however, as the spheroid of model 
3 is sufficiently compact to avoid widespread stripping. 

Similar to the case of the disk, the unbound spheroid 
particles fall into orbits on both the leading and trailing 
side of the orbit, with a sig nificant fraction being en- 
gulfed by the LMC (Fig- 



panels). By t 
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panel), the MW has greatly elongated the leading and 
trailing streams. The second close passage of the MCs 
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Fig. 12. — Close-up of the final SMC spheroid distribution pro- 
jected onto the sky for model 1 (top panel, cyan contours), model 
2 (middle, black contours), and model 3 (bottom, magenta con- 
tours). Contours are drawn at projected particle counts of 10 x 
where x =(0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5), and the filled regions 
range from light grey to dark grey as x increases. The colors of 
the contours merely indicate the identity of each model (i.e., the 
contours represent the same densities, regardless of color). The on- 
sky coordinates are Magellanic longitude Lms an d latitude Bms as 
defined by Nidever et al. (2008). 

at t = —0.26 Gyr successfully removes material from all 
three spheroid models, though model 3 once agai n su ffers 
the least amount of stripping. At t = Gyr (Fig [IT] bot- 
tom right panel), the MCs have arrived to their present- 
day positions and velocities, and the spheroid particles 
have reached their final distribution. As can be seen, the 
leading and trailing streams contain particles from mod- 
els 1 and 2 only, but particles from all three models are 
stripped into the region between the LMC and SMC (i.e., 
the Magellanic Bridge region) . Because the spheroid im- 
itates an old stellar population (e.g., Harris & Zaritsky 
2006), the stripped structures can be taken as possible 
stellar analogues of the MS, LA, and Bridge. 
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Fig. 13. — The on-sky radial density profiles of the SMC spheroid 
for models 1 (cyan), 2 (black), and 3 (magenta) at the present 
day t = 0. To construct the profiles for each model, the spheroid 
particles are binned according to on-sky radius (bin size 0.1°) and 
azimuthal angle (bin size 51°). We use the on-sky location of the 
SMC center of mass in each model as the origin r = 0, and we 
require at least one particle on average in each bin. The solid lines 
indicate the average density at each radius (i.e., annular ring in 
the sky) for each model, and the shaded region gives the standard 
deviation across the azimuthal bins. Observational data (red circles 
and error bars) is taken from the red giant survey of Nidever et al. 
(2011). In order to facilitate comparison, the density profiles of 
the models have been rescaled such that the density of model 1 
matches the Nidever et al. (2011) datapoint at r f» 2° (see text for 
details). 

4.2. SMC Outer Halo 

In their recent survey of the stellar periphery of the 
SMC, Nidever et al. (2011; Nil) discovered an extended 
halo component expanding to at least r ~ 11° from the 
SMC. Although it is observationally unclear whether this 
newfound component is bound or unbound, we can test if 
the tidal distortions of the SMC spheroid in our adopted 
model are able to reproduce this important new obser- 
vation. Fig [12] gives a close-up view of the SMC for our 
three spheroid models. Models 1 and 2 (top and middle 
panel, respectively) exhibit an extended halo population, 
but such a component is largely absent in model 3 (bot- 
tom panel), which remains dense and compact. A bit of 
tidal distortion is evident for model 3 toward the Bridge 
region (i.e., to the lower left of the panel), but this can- 
not be described as an extended, azimuthally symmetric 
population as discovered by Nil. 

Further to their discovery, Nil show a distinct "break 
point" in the radial density profile of red giants occur- 
ring r ~ 7° from the SMC. This break point (i.e., sharp 
change in slope) is indicative of tidal distortions (e.g., 
Munoz et al. 2008) and should be reproduced in tidal 
models of the SMC. In Fig[l3]we compare the Nil radial 
density data against the predictions of our three spheroid 
models, where lines give the average radial density of 
each model and shaded regions give the standard devia- 
tion of different azimuthal bins. In other words, Fig [l~3] 
provides not only the radial density but also its degree 
of symmetry in different directions of the sky. Nil ob- 
serve the SMC halo to be roughly symmetric in different 
azimuthal directions, and Fig [13] shows that this prop- 
erty is suitably reproduced (i.e., the standard deviations 
of each model are comparable to the spread in observed 
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data at each radius). 

In order to facilitate comparison with Nil, the models 
have been scaled to the observed data such that the den- 
sity of model 1 matches the Nil datapoint at r w 2°. 
The same scaling is applied to all three models. We 
stress that our simulation neglects stellar evolution and 
cannot provide unique predictions for the density of red 
giants as observed by Nil. The scaling of the density 
profiles in Fig [13] merely provides a visual comparison, 
and we justify this scaling as follows. We cannot esti- 
mate the density of stars without arbitrarily assuming 
mass to light ratios, initial mass function, star formation 
history, etc., which is not possible in the present frame- 
work. In addition, what is observed by Nil is only a frac- 
tion of the SMC stellar population, i.e., a portion of the 
red giant branch, which further separates the observed 
and simulated quantities. Lastly, there are magnitude 
limitations that should be considered for the observed 
stellar densities. As discussed in Nil, previous studies 
(e.g. De Propris et al. 2010) failed to detect the ex- 
tended SMC halo due to magnitude limitations, and it is 
feasible that future surveys utilizing deeper photometry 
than Nil may detect even more stars. Despite these dif- 
ficulties, the combined effect of these factors would likely 
change only the normalization of the radial density pro- 
file, so we choose to scale our models to the Nil data in 
Fig[l3]to facilitate c omp arison. 

Judging from Fig |13| models 1 (cyan) and 2 (black) 
successfully populatethe halo of the SMC at the large 
radial distances (up tor~ 11°) observed by Nil. In ad- 
dition, the models predict a "break population" (r > 4°) 
having a much shallower slope for the radial density in 
comparison to the slope of the inner SMC population 
(r < 4°). Moreover, the simulated slopes compare fa- 
vorably to the observed values. However, the densities 
for models 1 and 2 do not agree with the Nil data at 
all radii, and the predicted position of the break point 
is at a smaller radius than observed (r ~ 7°). Neverthe- 
less, models 1 and 2 are far more successful than model 
3 (magenta), which fails to reproduce an extended halo 
population for the SMC. Even though the slope of its 



radial profile in Fig 
in models 1 and 2, t 



13] abruptly changes at r ~ 4° as 
;here are an insignificant n um ber of 
particles stripped away in model 3. Indeed, Fig 12 shows 
that few particles are able to populate the outer regions 
of the SMC in model 3, and the distribution furthermore 
does not resemble the azimuthally symmetric halo dis- 
covered by Nil. 

Accordingly, we consider an extended spheroid like 
models 1 and 2 to better represent the tidal evolution 
of SMC. Although the compact spheroid of model 3 is 
unable to reproduce the observations of Nidever et al. 
(2011), there is still the potential for unexplored interme- 
diate cases 2.5 kpc < i? sp h,SMC <5.0 kpc. In this sense, 
a compact spheroid is disfavored but not ruled out. 

Before moving on, we compare our models to three 
more observational details of the SMC stellar halo: (i) 
Nil find that the population of red giants at radii 3° < 
r < 7.5° has a center that is offset by ~ 0.59° from the 
center of the inner stellar distribution; (ii) Nil estimate 
the radial scale length of the SMC to be 0.8° (1.0°) when 
using spherical (elliptical) model fits; (hi) The ellipticity 
of the SMC halo at 3° < r < 7.5° is estimated by Nil to 
be e ~ 0.1, which is less flattened than the inner parts of 



the SMC (e « 0.3 at r < 3°, Harris k Zaritsky 2006). 

Comparing against the offset of centers (i), we find off- 
sets of - 0.59° for model 1, ~ 0.33° for model 2, and 
~ 0.20° for model 3, where the offsets are calculated 
between the outer population (3° < r < 7.5°) and in- 
ner population (r < 3°) in each model. The cause of 
this offset can in general be attributed to the different 
responses to tidal perturbations between the inner and 
outer parts of the SMC, but perspective effects must be 
considered as stressed by Nil. To compare against the 
observed scale length (ii), we calculate the radial scale 
length in each of our spheroid models, which is roughly 
equal to 0.8Rh, where Rh is the half-mass radius at the 
present day. We find scale lengths of 0.6°, 0.4°, and 0.2° 
for models 1, 2, and 3, respectively. Considering the 
ellipticity (iii), we find that the models exhibit the op- 
posite behavior as observed, with greater spherical sym- 
metry at small radii (where the gravitational potential 
is strongest) and more highly flattened distributions at 
large radii (where the potential is wea ker and tidal dis- 



turbances are stronger). From Fig 121 one can estimate 
that all models are roughly spherical e « in the in- 
nermost region, and in the outer regions the ellipticity 
approaches e « 0.3 for models 1 and 2, and e < 0.1 for 
model 3. In summary, we find points (i) and (ii) reinforce 
our preference for the extended spheroid models 1 and 
2, whereas the discrepancies in ellipticity for point (iii) 
await the attention of future work. 

4.3. On-sky Distribution and Distances 

In Fig [14] we project the final distribution of the 
spheroid models onto the sky using a ZEA projection 
centered on the south galactic pole. Comparing against 
Fig [5] we see that the tip of the MS is offset from the 
stellar streams of models 1 and 2 (cyan and black, respec- 
tively). This result is not unnatural because material re- 
moved from a rotating disk (e.g., the MS) may evolve 
differently to the material removed from a pressure- 
supported spheroid (e.g., the stellar stream). Past ob- 
servational attempts to locate stars along the MS have 
largely failed, such as Recillas-Cruz (1982; "plus signs" 
in Fig_ M| and Guhathakurta & Reitzel (1998; "cross" in 
Fig [l4|, who could not find stars at the tip of the MS. 
Because the location of these survey fields are offset from 
the location of the streams of models 1 and 2 in Fig 14 we 



suggest that the Magellanic system may possess a stellar 
stream as long as it is displaced from the gaseous MS. 
However, the existence of the stream can be challenged 
by surveys such as Bruck & Hawkins (1983; rectangle in 
Fig [T4| , whose observations at the base of the MS overlap 
with the predicted streams of models 1 and 2. Interest- 
ingly, the possible offset between the stripped gas and 
stars is echoed in the results of Bruck & Hawkins (1983), 
because they observed a patch of stars at the base of the 
MS that did not correlate with the peak HI regions. 

If it exists, the stellar analogue of the MS may be diffi- 
cult to observe. As shown in the top panel of Fig [15] the 
distances to the stellar stream (Lms < —20°, right hand 
side of Fig 15 ) are quite large, with many stars beyond 
100 kpc at the tip of the stream. Even though the trail- 
ing stream of model 1 is quite dense (1.2 x 10 3 M /deg 2 , 
on average), the stream of model 2 is much more diffuse 
(1.8 x 10 2 Mo/deg 2 on average). Even though it is not 
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Fig. 14. — The final distribution of spheroid particles projected 
onto the sky for models 1 (cyan), 2 (black), and 3 (magenta). The 
location of survey fields that failed to find stars correlated with the 
Magellanic Stream are shown as plus signs (Recillas-Cruz 1982), 
a cross (Guhathakurta & Reitzel 1998) and rectangle (Bruck & 
Hawkins 1983). The coordinates are galactic longitude I (straight 
lines) and latitude b (concentric circles) shown in a ZEA projec- 
tion centered on the South Galactic Pole. The orbit trails of the 
SMC (green) and LMC (red) are also shown, as are their current 
locations (circles). 

clear what fraction of these simulated densities would 
correspond to observable stars, it is nevertheless certain 
that the stream of model 2 would be considerably more 
challenging to detect than that of model 1. The stellar 
analogue of the LA (Lms > 10°, left hand side of Fig 15 1 
may also be difficult to observe due to confusion with the 
disk of the MW around b w 0. 

In contrast, the stripped stars in the Magellanic Bridge 
area should be relatively straightforward to observe. 
Stars have indeed been located in the Bridge (e.g., Irwin 
et al. 1985), but they are young and blue, rather than old 
and red as would be expected from the spheroid popula- 
tion. Indeed, Harris (2007) has suggested that the mate- 
rial stripped into the Bridge was almost purely gaseous, 
and any stars currently in the Bridge would have formed 
in situ. This may pose a problem in the current frame- 
work, as even our most compact SMC spheroid (model 
3, i? S ph,SMC = 2.5 kpc) is easily stripped into the Bridge 



region, where the average density is 3.6 x 10 3 M Q /deg . 
The corresponding Bridge densities for models 1 and 2 
are 5.7 x 10 4 M Q /deg 2 and 2.6 x 10 4 M /deg 2 , respec- 
tively. We c omment further on this topic in the Discus- 
sion Sec 15.11 



4.4. Kinematics and SMC Velocity Profile 



The bottom panel of Fig 15 gives the radial velocity 
profile of the stripped spheroid material. The profile for 
the trailing stream of models 1 and 2 (Lms < —20°, right 
hand side of Fig 15 ) is largely coincident with that of the 
stripped disk particles (Fig[7j lower panel), and follows 
the kinematic trends of botn the MS and Horn. The 
trailing streams of models 1 and 2 also exhibit an ap- 
parent "pile- up" of particles with vlsr ~ —200 km s _1 
between Lms = —100° and —150°, which mirrors the ve- 
locity inflection at the tip of the gaseous MS (Nidever et 
al. 2010), but it is important to note that these particles 
are distributed over a large range of distances (Fig 15 
upper panel). 

The velocity profile for the leading stream (Lms > 10°, 
left hand side of Fig 15 1 is similar to that of the simu- 
lated LA (Fig[7j lower panel), but it is considerably more 
extended to lower velocities at Lms > 70°. The particles 
for model 3 are st ripp ed only into the Bridge region, as 
is clear from Fig 15 We remind the reader that the 
tidal evolution of the disk failed to explain the small ro- 
tational amplitude of the intermediate-age sta rs o f the 
SMC (i.e., as measured by HZ06; see Section 3.6). To 
reveal the velocity structure of stripped SMC stars (e.g., 
Fig 15 ) future models must therefore represent the SMC 
as anspheroid plus disk" as in the present work, depart- 
ing from the traditional "pure disk" models of the past 
(e.g., Connors et al. 2006; Besla et al. 2010). 
The SMC velocity profile for each of the spheroid mod- 



els is given in Fig 



10 



Models 1, 2, and 3 (panels b, c, 
libit a velocity gradient of ~10 km 



d, respectively) all ex 

s _1 deg -1 across the face of the SMC. This corresponds 
nicely to the observations of Harris & Zaritsky (2006), 
who report a similar gradient of 8.3 km s _1 deg -1 . In 
their analysis, HZ06 determine that the maximum pos- 
sible rotation velocity that could be derived from this 
velocity gradient is significantly less than the velocity 
dispersion. For this reason, HZ06 conclude that the the 
stellar component of the SMC is a dispersion supported 
spheroid. All three of our models are consistent with this 
analysis, but one may still ask why the velocity gradient 
is imprinted into the SMC. Harris & Zaritsky (2006) 
consider that that the observed velocity gradient may be 
caused by a projection effect of the SMC space velocity 
along the line of sight. That is, because the SMC oc- 
cupies a large area on the sky, the projection of radial 
velocity at different regions of the SMC will necessarily 
vary (e.g., van der Marel et al. 2002). 

In contrast, we consider that the velocity gradients in 
our simulations are due to tidal disturbances. As we 
have shown previously, much material is stripped away 
from the SMC during the close passages of the LMC at 
t = —1.97 Gyr and t = —0.26 Gyr. What we emphasize 
here is that these close encounters also strongly affect 
the particles that remain within the SMC. In particular, 
each close passage of the LMC effectively reorganizes the 
spheroid particles along a velocity gradient. If the parti- 
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Fig. 15. — The heliocentric distances (top panel) and t>LSR velocities (bottom panel) of the spheroid particles for models 1 (cyan), 2 
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cles are pulled beyond the tidal radius of the SMC, as in 
the case of models 1 and 2 for the t = —1.97 Gyr passage, 
then the velocity gradient determines how the leading 
and trailing streams separate from the SMC body. Even 
though the particles of model 3 remain bound to the SMC 
during this interaction, the entire spheroid is neverthe- 
less imprinted wit h a v elocity gradient. As shown at the 
present day in Fig 10 (i.e., at t = Gyr), the spheroids 
exhibit velocity gradients that are the signature of the 
most recent LMC passage at t = —0.26 Gyr. 

4.5. Accretion onto the LMC 

Fig [l6| shows the total mass falling into LMC (r < 7.5 
kpc) aifunction of time for the SMC disk (model 1) and 
the SMC spheroid (models 1, 2, and 3). We stress 
that Fig 16 merely gives the particles passing through 
the LMC at each time step, whereas the fate of each par- 
ticle (i.e. bound or unbound to the LMC) depends on 
its energy. The mutual close encounters of the MCs at 
t = —1.97 Gyr and t — —0.26 Gyr bring their centers of 
mass to within 6.0 kpc and 6.6 kpc, respectively, which 
explains the large temporary peaks in Fig |16[ However, 
a portion of the SMC particles are truly transferred to 
the LMC following each epoch of interaction. Up to the 
present day (t = Gyr), the total mass transfered to 
the LMC equals - 3.1 x 10 7 M from the SMC disk, 
~ 5 x 10 6 M from the extended spheroid (model 1), 
~ 2x 10 6 Mq from the intermediate spheroid (model 2), 



and ~ 2 x 10 5 M Q from the compact spheroid (model 
3). For each component, approximately half of the mass 
was captured by the LMC following the first interaction 
at t = —1.97 Gyr and the other half following the sec- 
ond interaction at t = —0.26. As can be seen in Fig 16 
the transfer of mass from the disk exhibits a sharp rise 
from - 10 6 M at t = -1.75 Gyr to - 1.3 x 10 7 M at 
t = —1.25 Gyr, in contrast to the comparatively flat dis- 
tributions of mass transfer from the spheroids. This can 
be attributed to the strong rotation of the disk particles 
as they fall gradually into the LMC with considerable 
angular momentum. It also explains the slight decrease 
of mass transfer from the disk between t — —1.25 Gyr 
and t — —0.8 Gyr because some particles that initially 
fall within 7.5 kpc of the LMC subsequently exceed this 
radius by virtue of their large angular momentum. 

Recently, Olsen et al. (2011) have discovered a popula- 
tion of AGB stars within the LMC that possess peculiar 
kinematics and metallicities. They suggest that these 
stars may have originated in the SMC and were accreted 
at some point during the recent dynamical interaction 
with the LMC. Our models provide a concrete realization 
of this scenario. In particular, because AGB stars would 
be represented by the old type stellar population of the 
SMC spheroid, the transfer of spheroid particles to the 
LMC reproduces the main result of Olsen et al. (2011). 
Fig [TT] gives the i^lsr velocity profile of the particles cap- 
turecTby the LMC (i.e., r < 7.5 kpc at t — Gyr) from 
the SMC disk (panel a), and the spheroids of models 1, 2, 
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Fig. 16. — The total mass falling within 7.5 kpc of the LMC as 
a function of time, shown for the SMC disk of model 1 (solid line), 
and the SMC spheroid of models 1 (long dash), 2 (dash), and 3 
(dotted). The two peaks correspond to the two close passages of 
the LMC and SMC. 

and 3 (panels b, c, and d, respectively). The properties 
of the accreted spheroid particles are broadly consistent 
with the results of Olsen et al. (2011), namely because 
the kinematics are distinct from that of the LMC disk. 
However, the exact kinematical signature of the Olsen et 
al. (2011) population, i.e., counter-rotation with respect 
to the LMC disk, is not reproduced. Considering that 
the population of stars that Olsen et al. (2011) attribute 
to the SMC comprise a small but significant percentage 
> 5% of their observed sample, the mass transfer sce- 
nario of model 3 (Fig 17 panel d) may be insufficient 
to explain the observation. This would seem to argue 
against a compact spheroid for the SMC. In Sec |5.1| we 
consider further implications of the Olsen et al. (2011) 
observation in the context of other results in the litera- 
ture. 

5. DISCUSSION 

5.1. The Possibility of Tidally Stripped SMC Stars 

Numerous surveys have shown that stars do not ex- 
ist in large numbers along the MS and furthermore 
that they do not correlate with peak HI regions (e.g., 
Recillas-Cruz 1982; Bruck & Hawkins 1983). Addition- 
ally, Guhathakurta & Reitzel (1998) use the null observa- 
tion of stars in their 5' x 7' field to derive an upper bound 
on the ratio of stellar to gaseous mass M^/Mm saO.l at 
the tip of the MS. Is the present model consistent with 
these constraints? We must be careful to consider two 
distinct populations when discussing the tidal evolution 
of SMC stars: (i) the stellar component of the disk, and 
(ii) the spheroidal component of old-type stars. As seen 
in Fig[l4j stars stripped from the spheroid do not coin- 
cide with the MS tip, nor is there a significant overlap 
with the probes of Recillas-Cruz (1982) or Guhathakurta 
& Reitzel (1998). In addition, there would be no reason 
for stars stripped from the spheroid to correlate with the 
high-density regions of the MS, w hich is consistent with 
the above observations (Sec 4.3 1. Stars stripped from 



the disk, however, would be interspersed among the gas 
clouds of the MS and would likely coincide with peak 
HI regions, which poses a problem if the stripped stars 
are too numerous. In other words, stars stripped from 



■j 



© 



1 1 1 1 I ' 1 1 .' I. .' .' 1 1 I 
■••' '•:'U 



® 



i f I- 



v 

1 . 

. ; 



5 
B MS (deg) 



5 0-5 
Bms (deg) 



185 



220 



255 
VLSR(km/s) 



290 



325 



Fig. 17. — The dlsr velocity field of particles accreted to within 
7.5 kpc of the LMC. The particles originate from (a) the SMC disk 
of model 1, and (b-d) the SMC spheroid of models 1-3, respectively. 
Each panel is centered on the LMC, and the dashed circle outlines 
a radius of 3 kpc. The on-sky coordinates are Magellanic longitude 
Lms an d latitude Bms as defined by Nidever et al. (2008). The 
velocity range has a lower bound of 185 km s — 1 (blue), an upper 
bound of 325 km s _1 (red), and is centered on the i>lsr velocity 
of the LMC (white). 



the SMC spheroid may be permitted by current observa- 
tional constraints, but those sourced from the SMC disk 
face stronger restrictions. 

To address this concern, previous models have assumed 
that the stellar component of the SMC disk was confined 
to the central regions where tidal stripping is inefficient 
(Yoshizawa & Noguchi 2003; Besla et al. 2010). In the 
present work, we naively associated t he col lisionless disk 
particles of our model with gas (Sec 2.3.2), but here we 
consider the possibility that some portion of the disk, 
particularly the cen tral regions, may pertain to stellar 
material. In Fig [18] we plot the number of disk particles 
as a function of their radius within the initial SMC disk 
(solid line), and we also show the fraction that eventu- 
ally constitute the MS (dashed) and the LA (dotted). 
The particles are selected to be in the MS or LA based 
on the criteria given in the caption of Fig [l8j Most of 
the stripped particles originate in the outer half of the 
SMC disk, with the average radius being 3.04 kpc and 
3.03 kpc for the MS and LA, respectively. If we assume 
that the disk mass within r = 1.75 kpc (2.0 kpc) is com- 
pletely stellar and beyond this radius is strictly gaseous, 
then the stellar to gas mass ratio is 0.09 (0.19), which is 
broadly consistent with the upper limit of ~0.1 measured 
by Guhathakurta & Reitzel (1998) at the MS tip. This 
would imply that the initial SMC gas disk (r =5.0 kpc) 
had an extent which was 2.5 to 2.85 times that of the 
stellar component in the disk, which is not uncommon 
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Fig. 18. — Number of particles plotted as a function of their orig- 
inal radius within the SMC disk. The solid line corresponds to all 
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(LA), respectively. The MS is selected as all particles satisfying 
Lms < — 25 that are not in the Horn, and the LA is selected as all 
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for dwarf galaxies (Swaters et al. 2002). 

Even though there is no direct evidence for a stellar 
counterpart to the MS, numerous studies indicate that 
the outer halos of the LMC and SMC are much more 
expansive than previously appreciated (e.g., Kunkel et 
al. 1997b; Noel & Gallart 2007; Saha et al. 2010; De 
Propris et al. 2010), extending more than 20° from the 
LMC (Muhoz et al. 2006; Majewski et al. 2009) and 
more than 10° from the SMC (Nidever et al. 2011). It 
would seem natural for these extended halos to be ac- 
companied by diffuse tidal streams of stars, and there 
are hints that such structures may indeed exist. For ex- 
ample, stars observed in the foreground of Carina may be 
the first evidence of tidal stripping from the SMC and/or 
LMC stellar halo (Muhoz et al. 2006). In addition, the 
radial density profiles of the LMC stellar periphery (Ma- 
jewski et al. 2009) and that of the SMC (Nidever et al. 
2011, hereafter Nil) both exhibit a sharp "break" be- 
yond which the slope of the density profile changes. Such 
a break is present in many other tidally evolved systems 
and may indicate the presence of tidal tails, as in the 
Sagittarius dwarf galaxy (Majewski et al. 2003), the Ca- 
rina dwarf spheroidal galaxy (Muhoz et al. 2006), and 
numerous globular clusters (Chun et al. 2010). However, 
the observation of a break population cannot by itself es- 
tablish the existence of tidal structures, because it may 
also be consistent with a bound stellar halo (e.g., Muhoz 
et al. 2008). 

In the present study, we are better able to reproduce 
the SMC break population (e.g., Fig 13 1 when we adopt 
an extended stellar spheroid (-R S ph,SMC 5 kpc) rather 
than an initially compact one (i? S ph,SMC = 2.5 kpc). 
In turn, the model predicts that diffuse tidal features 
are stripped away from the extended spheroid during its 
tidal evolution, forming stellar structures analogous to 
the gaseous MS, LA, and Bridge. In other words, we have 
used the Nil observations of the SMC stellar periphery 
to support the hypothesis that stellar tidal structures 
exist in the Magellanic system. Future wide-field obser- 
vational surveys such as the Southern Sky Survey by the 
SkyMapper telescope (Keller et al. 2007) would be ideal 



for detecting the stellar counterparts to the MS and LA, 
if they indeed exist, and establishing a connection to the 
outer halo of the SMC. 

The present prediction that stars should be stri pped 
into the Bridge region from all spheroid models (Sec 4.3 1 
is in conflict with Harris (2007), who observed the den- 
sity of old-type stars in the Bridge to be quite low. Young- 
blue stars have indeed been observed in the Bridge, but 
they were not tidally stripped from the SMC but rather 
were formed in situ following the epoch of tidal strip- 
ping (Irwin et al. 1985). It would appear that a much 
more compact spheroid (i? sp h,SMC ^ 2.5 kpc) would be 
needed to avoid stripping into the Bridge, but such a 
spheroid would unlikely be able to reproduce the ex- 
tended halo structure observed by Nil. Furthermore, 
the recent contention of Olsen et al. (2011) that the 
LMC has accreted a population of AGB stars from the 
SMC also creates tension with Harris (2007), because 
we find that tidal stripping through the Bridge region 
provi des a me chan ism of such a mass transfer (see sec 



tions 4.3 and 4.5). Another curiosity follows from the 
large line-of-sight depth of the SMC, which has been ex- 
plained in the present model as a pr ojection of the tidally 
extended Counter-Bridge (Sec 3.5). However, this ex- 
planation relies on stars being tidally stripped into the 
Counter-Bridge, which would logically imply that stars 
should be stripped into the Bridge as well. Regardless of 
these details, it is clear that many recent studies (includ- 
ing the present work) have placed a renewed interest in 
the tidal history of the SMC stellar population. 

5.2. Hydrodynamical Interactions with the Milky Way's 

Hot Halo 

The present work shows that gravitational interactions 
alone can explain the morphology and kinematics of the 
MS and Bridge, but the LA on the other hand is not 
well explained in the current framework. Here we con- 
sider the second-order effects of including hydrodynam- 
ics and its potential for reproducing the LA. The major 
hydrodynamical interactions to be considered are encom- 
passed by the collision of neutral gas clumps (e.g., within 
the MS and LA) with the tenuous hot halo of the Milky 
Way. The drag induced by this interaction is termed 
"ram pressure" and is proportional to pv 2 where p is 
the density of the hot halo at the location of the HI gas 
clump, and v is the relative velocity between the hot 
halo and the clump (e.g., Gunn & Gott 1972; Westmeier 
et al. 2011). As a simple approximation, the density p 
can be considered a function of galactocentric radius r, 



with reasonable choices scaling as p 



(isothermal) 



or p ~ r~ 3 (NFW), and the velocity v can be taken as 
the galactocentric velocity if we assume the hot halo is 
non-rotating. These simple assumptions will allow us to 
calculate the relative effect of ram pressure on the MS 
and LA throughout their orbit about the MW. 

Fig 19 gives the galactocentric distance (top panel) 
and average velocity (bottom panel) throughout the en- 
tire orbital evolution of the present model. The MS and 
LA are selected as in Fig [l8j but the LA is split into 
two branches (LAI and LA2) according to their distinct 
evolutionary tracks. As can be seen, the orbits of the 
LA branches are generally faster and at smaller galactic 
radii than that of the MS, and the corresponding effect 
of ram pressure will accordingly be larger. For example, 
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at t = —0.3 Gyr, one can calculate that the relative effect 
of ram pressure on the LA should be at least an order of 
magnitude larger than for the MS. Guided by this sim- 
ple calculation, we consider hydrodynamical effects to be 
relatively minor for the MS but very important for deter- 
mining the morphology and kinematics of the LA. This 
contention will be explored with future, more sophisti- 
cated N-body models. In our previous work (DBllb) 
we argued that the distances and velocities of the LA 
can decrease in response to a simplistic prescription for 
ram pressure, and that such effects would produce an 
improved agreement with observations (e.g., McClure- 
Griffiths et al. 2008). Under these expectations, strong 
ram pressure would likely cause the LA to sink to lower 
radii subsequent to t = —0.3 Gyr in the present model. 
However, because we adopt collisionless dynamics in the 
present work, the LA is instead predicted to ris e to much 
larger radii subsequent to t = —0.3 Gyr (Fig 
this reason, we cannot consider the present col 
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model to adequately capture the evolution of the LA. 
This explains, at least partially, why the predicted on- 
sky location in Fig[5]is quite different to that observed. 

Hydrodynamical effects at the interface of the hot halo 
are also responsible for the shredded, ionized filaments 
of the MS (Stanimirovic 2008; Westmeier & Koribal- 
ski 2008), the disjointed arrangement of LA clouds and 
their head-tail structure (Bruns et al. 2005), and the 
anomalous Ha emission within the MS (Putman et al. 
2003b; Bland-Hawthorn et al. 2007). All of these in- 
teractions lead to a decrease in the overall HI mass by 
transforming the neutral gas to an ionized component. 
The observed HI masses of the MS and LA are 1.2 x 10 8 
Mq and 0.3 x 10 8 M (Bruns et al. 2005), respectively, 
and the corresponding predicted masses of the present 
model are 1.2 x 10 8 M Q and 0.2 x 10 8 M , where we 
have used the definitions of Fig [18] to select the MS and 
LA. Note that Putman et al. (2003a) quote a larger MS 
mass of 1.9 x 10 8 M due to a larger region of selection 
for the MS. Even though the predicted masses compare 
well with the measured values, the inclusion of hydro- 
dynamical interactions would have the unwanted effect 
of reducing the predicted HI masses. This discrepancy 
could be combated by increasing the initial SMC mass, 
but this would alter the orbital history and may not be 
favorable for the formation of the MS. Another possi- 
ble solution is to change the mass ratio of the gas disk 
while keeping the total SMC mass fixed. We leave this 
nontrivial issue for future work. 

5.3. Mass Transfer from the SMC to the LMC 

The present models have clearly shown that gas and 
stars in the SMC can be transferred to the LMC via 
tidal interaction over the last ^2 Gyr. In particular, 
two epochs of significant mass transfer are predicted, one 
following the first strong interaction between the MCs 
at t = —1.97 Gyr, and the other following the second 
interaction at t = —0.26 Gyr. Here we consider how this 
predicted scenario compares with observational evidence. 

5.3.1. Gas Transfer to the LMC 

Through collisionless dynamics alone, our model pre- 
dicts that the transfer of gas to the LMC (i.e., transfer 
from the SMC disk) totaled - 1.3 x 10 7 M Q by t = -1.25 
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Fig. 19. — The time evolution of the average radius (top panel) 
and velocity (bottom panel) of the SMC (green solid line) and 
its various stripped components, including the Magellanic Stream 
(MS; blue long dashed line) and the two branches of the Leading 
Arm (LAI and LA2; red dashed and r ed d otted line, respectively). 
The MS and LA are selected as in Fig |18| with an additional divi- 
sion of the LA according to heliocentric distances and ulsr veloc- 
ities. LAI is the branch having larger distances and slightly larger 
radial velocities; it is colored yellow-orange in Fig [7] and originates 
in tidal arm A. LA2 has smaller distances and lower radial veloc- 
ities; it is colored black in Fig [jj and originates in tidal arm B. 



Gyr and ~ 3.1 x 10 7 M by the present day. However, 
these figures would increase, perhaps dramatically, once 
hydrodynamics is considered, because collisional forces 
would precipitate the infall of much material that other- 
wise passes straight through the LMC at t « —2.0 Gyr 
and t ss —0.25 Gyr (i.e., the two temporary peaks in 
Owing to the fact that the gas from the SMC 



Fig [16|. 

is relatively metal-poor, the accretion of this gas onto 
the LMC would provide fuel for the formation of low- 
metallicity stars (Bekki & Chiba 2007). The prediction 
that such a process has occurred is supported by two sep- 
arate pieces of observational evidence: (i) the chemical 
enrichment history of the LMC, and (ii) the unusually 
low abundance of nitrogen in young stellar populations 
of the LMC. 

we used the correlated star for- 



As described in Sec 2.1 



mation histories of the LMC and SMC (Harris & Zaritsky 
2009) to constrain the timing of their close tidal interac- 
tions for the past ~ 2 Gyr. But how did the metallicity 
of the LMC evolve during this time? The chemical en- 
richment history of the LMC is provided by numerous 
authors (e.g., Figure 20 of Harris & Zaritsky (2009); Fig 
1 of van Loon et al. 2005; and the ongoing VMC sur- 
vey, Rubele et al. 2012), which shows that the average 
metallicity of the LMC has generally increased over its 
~ 10 10 Gyr lifetime but has decreased during the most 
recent ~ 2 Gyr. Considering Fig 20 of Harris & Zaritsky 
(2009), the metallicity has decreased during two epochs: 
between ~ 2 to ~ 1 Gyr ago, and between ~ 500 Myr 
ago to the present day. These curious trends could possi- 
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bly be explained by widespread star formation from rela- 
tively metal-poor gas ([Fe/H]< —1), and the accretion of 
such gas onto the LMC is indeed predicted by the present 
model. This gas would likely admix into the LMC and 
effectively dilute the interstellar medium from which new 
stars would form. Notice also that the periods of decreas- 
ing metallicity in the LMC are remarkably well aligned 
with the two epochs of correlated star formation (~ 2 
Gyr ago and ~ 500 Myr ago; Harris & Zaritsky 2009). 
In fact, one can suggest that the decrease in metallicity 
appears to be a response to the bursts of star formation. 
Following the predictions of the present model, the two 
strong tidal encounters between the MCs (at t = —1.97 
and t = —0.26 Gyr) would have caused increased star 
formation in each galaxy, which explains the correlated 
bursts of star formation. And because the metal-poor 
gas from the SMC is transfered to the LMC in response 
to each tidal encounter, the stars in the LMC would sub- 
sequently form with relatively lower metallicities, which 
explains the decreasing trends in the chemical enrichment 
history. However, we should point out that the dominant 
mechanism for triggering the star formation in the LMC 
is not at all clear; it was likely a competition between the 
tidal interaction with the SMC and the gas infall from 
the SMC. 

Recent observations have confirmed that the chemical 
abundances of star-forming HII regions and very young 
stellar populations in the LMC have [N/H] ratios that 
are a factor of 6 — 7 lower than the solar value (e.g., 
Korn et al. 2002; Hill 2004; van Loon et al. 2010). On 
the other hand, other elements such as O and Ne are 
underabundant by a factor of only ~ 2 in these systems 
(Hill 2004). In other words, the young stellar popula- 
tions of the LMC exhibit curiously low abundances of 
nitrogen. To explain this anomaly, Bekki & Tsujimoto 
(2010) have invoked recent gas accretion onto the LMC, 
and they find that the most likely origin of this gas was 
from the SMC rather than, for example, high-velocity 
clouds (HVCs) in the Local Group. In particular, this 
gas would have originated in the outer parts of the SMC 
where metallicities were low owing to a negative metal- 
licity gradient (Bckki & Tsujimoto 2010). The present 
model provides a concrete realization of this scenario as 
long as we assume that the stripped gas is deficient in 
nitrogen. This assumption can be justified by the fol- 
lowing argument. Much of the gas that is accreted onto 
the LMC in our model has a similar origin as the MS, 
namely, as the tidally stripped material in the outer re- 
gions of the SMC disk, and it stands to reason that their 
chemical abundances should also be similar. Given that 
Fox et al. (2010) have observed very low nitrogen abun- 
dances in the MS, it is not unreasonable to assume that 
the gas transfer to the LMC was similarly nitrogen-poor. 
Under this assumption, the present model strengthens 
the conclusions of Bckki & Tsujimoto (2010), although 
it is unclear what amount of nitrogen-poor gas must be 
transfered to explain the properties of the young stel- 
lar populations of the LMC. We leave this quantitative 
investigation to our future work. 

5.3.2. Kinematically Peculiar Stars in the LMC 

Olsen et al. (2011) have recently investigated the 
distribution and kinematics of AGB stars in the LMC 
and found that about 5% of the stars in their sample 



have line-of-sight velocities that oppose the sense of ro- 
tation of the LMC disk. It is unlikely that these stars 
formed in the LMC, because in addition to having pe- 
culiar kinematics, they have much lower metallicities 
([Fe/H]= -1.25 ± 0.13) than the field stars of the LMC 
([Fe/H]= -0.56 ±0.02). Olsen et al. (2011) conclude 
that the most likely origin for these stars was the SMC, 
principally because De Propris et al. (2010) measured a 
similar metallicity of [Fc/H]= —1.35 ± 0.10 for the red 
giants in the periphery of the SMC. Two HI arms within 
the LMC are kinematically linked to the anomalous AGB 
population, and though these gaseous structures were 
previously identified as outflows (Staveley-Smith et al. 
2003), they are hypothesized by Olsen et al. (2011) 
to be infalling material. The present model reinforces 
this general scenario and illustrates a dynamical process 
whereby old-type stars in the SMC spheroid are tidally 
stripped and engulfed by the LMC. The model also pre- 
dicts that the stripped components of the SMC (both 
disk and spheroid) can exhibit rotation in the outer part 
of the LMC, possibly explaining the observed out-of- 
plane stellar polar ring (Kunkel et al. 1997a) and the 
counter-rotating stellar component in the LMC (Subra- 
manian & Subramaniam 2009). 

The transfer of stars from the SMC to the LMC occurs 
during two different epochs (t ss —2.0 Gyr and t w —0.25 
Gyr) following the close encounters between the MCs. 
The second epoch corresponds to the formation of the 
Magellanic Bridge, and the stripped stars from the SMC 
spheroid should indeed coincide spatially with the gas 
in the Bridge. Even the most compact spheroid model 
that we adopt is easily stripped into the Bridge region 
(Sec 4.3), but this prediction may conflict with previ- 
ous failed attempts to detect old-type stars within the 
Bridge (e.g., Harris 2007). We should point out that 
tidally stripping stars through the Bridge region provides 
an obvious mechanism of mass transfer to the LMC, but 
it is not the only mechanism in our model, since consid- 
erable mass is also transfered during the first interaction 
between the MCs at t « -2.0 Gyr. 

In the present model, the total stellar mass that the 
LMC accretes is at most ~ 5x 10 6 M Q , which comes from 
adopting an extended spheroid (Model 1, i? sp h = 7.5 
kpc). Though Olsen et al. (2011) claim an SMC origin 
for 5% of their sample, it is not immediately clear if the 
total mass of this population within the LMC should 
be estimated as 5% of the current LMC stellar mass 
- 3 x 1O 9 M (van der Marel et al. 2002). If we do make 
this simple assumption, then the total mass of SMC stars 
transfered to the LMC should be - 1.5 x 10 8 M Q . This 
figure is a factor of ~ 30 larger than predicted, and is in 
fact comparable to the total initial mass that we assume 
for the SMC spheroid, M sph = 2.7 x 10 s M Q . This may 
imply that the observational estimate of the mass frac- 
tion of accreted SMC stars is overestimated, or it may 
suggest that a different accretion/merger event for the 
LMC (e.g., minor merging of small satellite galaxies) is 
needed to explain the origin of the kinematically peculiar 
AGB stars. A final possibility is that the SMC may have 
had a significantly larger total mass before it began to in- 
teract with the LMC and Milky Way, which may be rea- 
sonable given that a recent study of the SMC's rotation 
curve proposed a larger mass of more than 6.5 x 10 9 Af Q 
(Bekki & Stanimirovic 2009). 
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We accordingly suggest that an estimate of the to- 
tal mass of the accreted SMC populations may help to 
constrain the original SMC mass prior to its interac- 
tions with the LMC. Although the only observations of 
possible accreted SMC populations is the AGB dataset 
of Olsen et al. (2011), ongoing photometric and spec- 
troscopic studies of planetary nebulae in the LMC will 
enable the determination of other accreted populations 
from the SMC. Moreover, these future studies will al- 
low us to discuss how the total masses of the accreted 
SMC populations might depend on stellar ages and abun- 
dances. In addition, the kinematically peculiar HI gas in 
the LMC can be used to give further constraints on the 
mass transfer process between the LMC and the SMC, as 
suggested by Olsen et al. (2011). In the future, we plan 
to run chemodynamical simulations to discuss how the 
gas and stars of the SMC can be transferred to the LMC 
and consequently influence its structure and dynamics. 

5.4. Diagnosis for recent MS formation models: Success 
and Failure 

The formation of the MS has been described by numer- 
ous dynamical models in recent years, and here we will 
point out (i) the key differences among their methodolo- 
gies, and (ii) their relative success and failure to explain 
the observed properties of the MCs, MS, Bridge, and 
LA. We divide the models into four groups: the tradi- 
tional tidal models (e.g., C06), the ram pressure strip- 
ping model of Mastropietro et al. (2005, hereafter M05), 
the "first passage" model of Besla et al. (2010, hereafter 
B10), and the present tidal model. Although there are a 
number of MS models based on test particle simulations 
(e.g., Ruzicka et al. 2010; DBlla), we focus here on 
the N-body and hydrodynamical simulations. Our com- 
parison is summarized graphically in Table [3] It should 
be stressed that we are simply attempting to highlight 
the relative strengths of the models rather than trying 
to discern which one is the "best." After all, each model 
will certainly be improved in the future as more realis- 
tic physical prescriptions are added and higher resolution 
simulations become accessible. 

Reproducing the orbital constraints of the MCs is one 
of the fundamentally important tasks for any model of 
the Magellanic system. The present work adopts the 
values given in Table [T] for the on-sky positions, dis- 
tances, line-of-sight velocities, and proper motions of the 
MCs. Of these parameters, only the proper motions are 
poorly constrained, although recent measurements have 
suggested that the low-velocity LMC orbits adopted by 
C06 and M05 may not be realistic. The high-velcoity 
LMC orbit adopted by B10 is justified by the proper 
motion measurements of Kallivayalil et al. (2006a), and 
the LMC orbit of the present model is supported by 
Vieira et al. (2010). While all models reproduce the 
LMC orbit (barring the difference in proper motions), 
the same is not true for the SMC. Perhaps for simplic- 
ity, M05 neglect the SMC entirely in their model, focus- 
ing only on the interaction between the LMC and the 
Milky Way. Because the present model shares the same 
methodology as C06, namely the backward orbit inte- 
gration scheme, the SMC orbit can be easily reproduced 
by suitable choice of parameters. The task is nontrivial, 
however, for the numerical method adopted by B10, and 
consequently their orbit for the SMC has the incorrect 
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The symbol "— " means that the model does not explain or does 
not address the given property; " v / " means that the model explains 
the property; and "(v0" means that the model can explain the 
property onl y som ewhat. Choices in each category are clarified in 
the text, Sec |5.4| 

1 "C06" is representative of the traditional tidal models: Connors 
et al. (2006); Yoshizawa & Noguchi (2003); Gardiner & Noguchi 
(1996). These models are very similar, as they adopt the same orbits 
for the MCs and use the same isothermal halo for the MW. "M05" 
stands for Mastropietro et al. (2005). "B10" stands for Besla et al. 
(2010) as well as Besla et al. (2012). 



on-sky position, line-of-sight velocity, and proper motion 
(see Besla et al. 2012). 

The underlying reason that B10 fail to reproduce the 
correct SMC orbit is their adoption of self-consistent in- 
teractions between the MCs, i.e., they represent both 
the LMC and SMC as ensembles of N-body particles. 
This choice certainly makes the B10 model more sophisti- 
cated than either C06 or the present work, which instead 
adopt a fixed potential for the LMC, but it also implic- 
itly introduces dynamical friction between the MCs, ren- 
dering pre-calculation of the orbits impossible. To clar- 
ify: in C06 and in the present work, present-day orbital 
constraints can be chosen as initial conditions for or- 
bit integrations, and the N-body evolution is guaranteed 
to follow these pre-calculated orbits; in self-consistent 
models, however, the N-body evolution will not obey a 
pre-calculated orbit. Instead, one must run many full- 
resolution models to test which initial conditions, if any, 
reproduce the present-day orbital constraints. Add to 
this the requirement that the orbit deliver a good tidal 
stripping model for the MS, and it is understandable why 
B10 were unable to satisfy all constraints. 

The B10 model is not fully self-consistent, however, 
because the Milky Way is modeled by a fixed potential. 
M05 self-consistently model the interaction between the 
LMC and Milky Way without resorting to external po- 
tentials, but we emphasize that the SMC does not exist 
in their simulation. The M05 model is the only one to 
include hydrodynamical interactions with the hot halo 
of the Milky Way, and the incorporation of this inter- 
actio n in to current tidal models is indeed a salient need 
(Sec 5.2). The ultimate goal to which all these models 



strive is the construction of a fully self-consistent model 
of the gravitational and hydrodynamical interaction be- 
tween the LMC, SMC, and Milky Way that (i) repro- 
duces all MC orbital constraints, (ii) models the Milky 
Way realistically, and (iii) explains the MS, Bridge, and 
LA. Even constructing a model that satisfies (i) and (ii) 
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has yet to be achieved, which means that we must use 
far simpler models (e.g., the current work) to investigate 
the formation of the MS. To be fair, the current model 
is not a satisfying representation of the Magellanic sys- 
tem, owing to various simplifying assumptions including 
the neglect of dynamical friction between the MCs and 
failing to represent the LMC as a self-consistent N-body 
system. 

C06 model the Milky Way as an isothermal sphere 
whereas the other three models in the current compar- 
ison choose NFW halos for the Milky Way. Owing to 
the prevalence of CDM cosmology, the C06 choice can- 
not be considered realistic. Unlike the present work in 
which a bulge and disk are also ascribed to the Milky 
Way, the potential adopted by BIO is a lone NFW halo 
with a virial mass of 1.5 x 10 12 M Q and a virial radius of 
240 kpc. These parameters fall in between those of the 
presently adopted NFW halo, M vir = 1.30 x 10 12 M , 
i? V ir = 175 kpc, and the halo of the "alternate" model, 
M vil = 1.90 x 10 12 M Q , i? vir = 269 kpc (see Table [TJ. 
These choices for the Milky Way should be compared 
with mass estimates of the Milky Way at large radii. For 
example, Gnedin et al. (2010) estimate that the total 
mass of the Milky Way within 80 kpc is Q.9~tH x 10 11 
M Q . To compare, the total Milky Way mass within 80 
kpc in our adopted model is 8.7 x 10 11 M and in our 
alternate model is 9.1 x 10 11 M . Both of these figures 
fall within the la uncertainty of the measured quantity. 
In this sense, the adopted potential of the Milky Way in 
the present work is realistic, and similar arguments hold 
for B10 and M05. 

The observed bifurcation of the MS into two distinct, 
parallel filaments is strongly reproduced in the present 
work, and we have put much effort into revealing the 
origin of this feature (Sec 3.3). C06 previously claimed 



to reproduce the MS bifurcation, but it must be noted 
that the C06 on-sky MS distribution does not convinc- 
ingly reproduce the location or morphology of the MS 
filaments. Neither do they identify a plausible dynami- 
cal origin for the bifurcation. Nevertheless, we can con- 
servatively state that the C06 model exhibits filamentary 
substructure within the MS, and moreover, that the M05 
and B10 models not not exhibit any such substructure. 
The MS in the B10 model exhibits a brief on-sky split, 
but the location and kinematics of this structure do not 
correspond well with observations. For this reason, we 
speculate that a filamentary MS may be attributable to 
tidal interaction over multiple passages of the Milky Way, 
as this is a shared property of the present work and C06 
but neither M05 or B10. 

Alternatively, Nidever et al. (2008) suggest that the 
filaments of the MS may have arisen from an internal 
"blowout" mechanism that expelled gas from the MCs. 
Using an intricate analysis of HI radial velocity profiles, 
they trace one of the filaments to an origin within the 
SMC, and they trace the other filament (and portions of 
the LA) to an origin within the LMC. The present model 
conflicts with these results, particularly because we are 
able to explain the dual filaments by invoking a tidal 
origin within the SMC only. Nevertheless, one of the 
MS filaments "appears" to pass through the LMC in our 
model as judged from the right hand panel of Figure [6j 
and we show in Sec 13.31 that this filament is connected 



to a branch of the LA which also "appears" to overlap 
with the LMC in Fig [6} Though superficially consistent 
with Nidever et al. (2008) in the sense that the MS fil- 
ament and LA seem to be connected to the LMC, this 
coincidence is merely a projection effect in our model. 
Regardless, we will consider possible hydrodynamical in- 
teractions between the MS/LA and the gas disk of the 
LMC in future work. 

A possible observational test of the "blowout" hypoth- 
esis of Nidever et al. (2008) is the measurement of metal- 
licity along each of the MS filaments and the LA. Because 
tidal models such as the present work assume the MS and 
LA originate from the SMC, the metallicities of these 
structures should be correspondingly low, as verified at 
the MS tip by Fox et al. (2010). In contrast, the blowout 
model requires much of the stripped gas to be enriched, 
as would be expected from gas originating in the LMC. 
In particular, the two filaments of the MS should have 
markedly different metallicities if one originated in the 
LMC and the other in the SMC as suggested by Nide- 
ver et al. (2008). Such metallicity measurements have 
not been conducted but would provide valuable insight 
into the origin of the HI structures of the Magellanic sys- 
tem. Another property of the MS filaments that bears 
discussion is that their kinematics exhibit sinusoidal os- 
cillations in the analysis of Nidever et al. (2008). This 
is not reproduced in the current model (Fig [7]) , and at- 
tempting to resolve this discrepancy will be a topic of 
future work incorporating hydrodynamical interactions. 

The density of the MS is observed to steadily decline 
across its length until it reaches low column densities at 
its tip (e.g., Putman et al. 2003a; Nidever et al. 2010). 
The traditional tidal models such as C06 are inconsis- 
tent with this gradient because they predict an increasing 
trend with a dense plume at the MS tip. Such a dense 
region exists in the present work as well, but it is much 
more localized and located further along the MS beyond 
its classical tip (see Appendix). Moreover, the present 
model does indeed predict a generally decreasing trend 
in density until the tip is reached. Similarly, the B10 
model predicts the correct column densities for the MS 
though not the correct gradient. Only the M05 model 
strongly reproduces the MS density gradient, implying 
that hydrodynamical interactions are needed to explain 
this property. DBllb showed that the dense plume at 
the MS tip in the traditional tidal scenario can effec- 
tively be removed if ram pressure from the Milky Way's 
hot halo is taken into account. Accordingly, it seems that 
the density of the MS cannot be explained by tidal mod- 
els without external gas dynamical interactions with the 
Milky Way. 

The recent observational results of Fox et al. (2010) 
indicate that the chemical abundances of the MS tip 

are rather low (e.g., [O/H] 1.00 and [N/H]< -0.44), 

much lower than the metal-rich gas of the LMC and even 
lower than the chemical abundances of the SMC (e.g., 
[O/H] ~ —0.66). Considering this measurement, Fox et 
al. (2010) suggest that the MS likely originated in the 
SMC rather than the LMC. This is a potentially serious 
problem for M05 because they assume an LMC origin for 
the MS, in contrast to the tidal models (C06, B10, and 
the present work) which assume an SMC origin. In these 
tidal models, the outer parts of the SMC are preferen- 
tially stripped to become the MS (e.g., Fig 18). If the 
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SMC has a negative radial gradient in metallicity (i.e., 
if the outer gas disk is metal-poor in comparison to the 
central body of the SMC), then each of the tidal models 
would be consistent with the chemical abundance mea- 
surements of Fox et al. (2010). 

The Magellanic Bridge is well reproduced in C06 and in 
the present work, but the gaseous structure which con- 
nects the MCs in the BIO model is displaced from the 
observed location of the Bridge. Because M05 exclude 
the SMC from their model, they cannot address the for- 
mation of the Bridge. The M05 model also fails to predict 
any stripped material on the leading side of the orbit, in 
contrast to the leading structures predicted by the tidal 
models. Indeed, the observed presence of the LA has 
long been considered to support a tidal origin for the MS 
rather than a ram pressure origin (Putman et al. 1998). 
However, the observed location, kinematics, and mor- 
phology of the LA have not been successfully explained 
by any dynamical model. This would seem to suggest 
that the tidal models need to incorporate hydrodynam- 
ical effects in order to acco unt for the properties of the 
LA, as discussed in Sec |5.2| 

The predicted properties of the leading material may 
be able to distinguish between the different tidal mod- 
els, even if the correspondence to observation is poor. 
The C06 model does the best job in capturing the on- 
sky "kink" in the LA as it extends away from the MCs, 
but the observed LA location is more than 30° away from 
the C06 prediction. This "kink" is not exhibited in the 
present work, though the tip of the predicted leading 
material does coincide with a portion of the observed 
LA (Figp5|. The situation is worse for the BIO model, 
as it predicts the LA to emanate from the far side of 
the LMC where no HI structures are observed. In ad- 
dition, the BIO model is unable to predict an elongated 
leading structure possibly because the potential of the 
Milky Way is unable to play a vital role within the short 
timescale (< 0.2 Gyr) of a first passage orbit. The cur- 
rent model predicts two tidal branches within the lead- 
ing material, in contrast to the single leading arms of 
BIO and C06. These multiple branches may correspond 
to the proliferation of cloudlets observed in the LA (e.g., 
Diaz & Bekki 2011a), but more extensive modeling is 
needed to seriously consider this point. 

One of the principal disadvantages of the C06 and B10 
models is that they cannot reproduce the observed stel- 
lar kinematics of the SMC. Neither can the M05 model, 
which excludes the SMC entirely. C06 and B10 represent 
the SMC as a rotating disk galaxy whereas the observa- 
tions of HZ06 clearly show that the red giants of the 
SMC are distributed in a dispersion-supported spheroid. 
In the present work, we have self-consistently repro- 
duced the kinematical differences between the gaseous 
and stellar components of the SMC by adopting a multi- 
component system composed of a disk (for the gas) and 
a spheroid (for the old-type stars). As shown previously, 
the tidal evolution of the spheroid can also explain the 
observed stellar structure in the periphery of the SMC 
and therefore be more consistent with observations. This 
suggests that future tidal models will need to adopt disk- 
plus-spheroid systems for the SMC rather than the tra- 
ditional choice of a pure disk. 



5.5. The Magellanic past before the formation of the MS 

In the present study we have proposed that the MCs 
have become binary companions only recently, and that 
their first strong tidal interaction ~ 2 Gyr ago was sig- 
naled by the formation of the MS. There are two general 
scenarios that can accommodate the MCs commencing 
their tidal interaction only recently. First, the MCs may 
have orbited as independent satellites of the Milky Way 
at early times but have only recently joined as a binary 
pair. This scenario is invoked by a number of studies to 
explain the formation of the MS (DBlla and the present 
work) and to solve the "Age Gap" problem for globular 
clusters of the LMC (Bekki et al. 2004). The second 
possibility is that the MCs were accreted onto the Milky 
Way as a binary pair, but the mutual separation between 
the MCs has become small only recently owing to the 
gradual decay of orbits from dynamical friction. One re- 
alization of such an orbit has been investigated by Besla 
et al. (2010; see also Besla et al. 2012). It is currently 
unclear which of the above two possibilities is more plau- 
sible for the origin of the MCs, due in part to the lack of 
fully self-consistent simulations (i.e., representing both 
the MCs and Milky Way with N-body particles) that 
successfully reproduce the MC orbital constraints (e.g., 
see Table 

Another important aspect of the adopted orbital model 
is that the MCs interact with the Milky Way during two 
pericenters, ~2.4 Gyr ago and the present day. We spec- 
ulate that the MCs must have executed at least two pas- 
sages about the Milky Way in order to reproduce the HI 
observations of the Magellanic system, particularly the 
MS bifurcation and the elongated LA (in contrast, con- 
sider the first passage models of Besla et al. 2010 and 
Besla et al. 2012). Recent cosmological simulations of 
the LMC orbital history indicate that a first passage or- 
bit is not necessarily the preferred scenario as long as the 
virial mass of the Milky Way is > 1.8 x 10 12 M© (Sales 
et al. 2011). In the present work, the "alternate" model 
for the Milky Way satisfies this requirement, as it has a 
virial mass of 1.9x 10 12 M (see Table!!]). The total mass 
of the "adopted" Milky Way model within r < 300 kpc is 
comparably large, totaling 1.7 x 10 12 M Q . Accordingly, 
we interpret our orbital model (and the preference for 
two pericentric passages about the Milky Way) as being 
broadly consistent with the cosmological simulations of 
Sales et al. (2011). 

To consider the orbital history of the MCs at early 
times, one should not use fixed potentials as in the 
present work because the orbit integrations are realistic 
only for the past ~ 4 Gyr or so. Recently, a number of 
studies have explored more reliable methods for address- 
ing the full orbital history of the LMC: Bekki (2011b), 
who uses a self-consistent treatment to investigate the 
accretion of the LMC onto the Milky Way from out- 
side its virial radius, as well as Sales et al. (2011) and 
Boylan-Kolchin et al. (2011, hereafter BK11), who iden- 
tify "LMC-type" galaxies within large-scale cosmological 
simulations (see also Lux et al. 2010). Bekki (2011b) 
determine that the orbital models adopted in the tradi- 
tional formation scenarios for the MS (e.g., GN96) are 
not inconsistent with the conclusions of BK11, and we 
arrive to a similar conclusion for the present model. The 
LMC has a number of orbital properties in the present 
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work including (i) it is strongly bound to the Milky Way, 
(ii) the orbital eccentricity is e ~ 0.5, and (hi) the LMC 
has been orbiting the Milky Way for at least 4 Gyr. Prop- 
erty (hi) implies that the LMC needs to be accreted onto 
the Milky Way more than ~4 Gyr ago, and BK11 deter- 
mine that such galaxies (i.e., LMC-analogs that cross the 
Galactic virial radius at tf c > 4 Gyr ago) may have small 
eccentricities and are more likely to be strongly bound. 
In particular, BK11 show that about 50% of such galax- 
ies also have eccentricities e smaller than 0.6. This result 
is consistent with the above properties (i) and (ii) of the 
present work. We therefore conclude that the present 
orbital model can be consistent with the predictions of 
cosmological simulations as reported by BK11. 

6. CONCLUSIONS 

We have constructed new N-body models for the tidal 
evolution of the SMC in a ~ 3 Gyr interaction history 
with the LMC and Milky Way, and we have thereby in- 
vestigated the formation processes of the MS, LA, and 
Bridge. The framework of the present investigation has 
a number of advantages over previous tidal models of the 
SMC (e.g., C06, BIO), including (i) our adopted orbital 
model was selected after conducting a large parameter 
space search centered on the measured proper motion 
values of Vieira et al. (2010), (ii) the Milky Way is rep- 
resented by a realistic potential having a bulge, disk, and 
NFW halo, and (hi) we investigated the stellar kinemat- 
ics and outer stellar structure of the SMC by adopting 
"disk plus spheroid" models in good agreement with ob- 
servations (e.g., Harris & Zartisky 2006; Nidever et al. 
2011). Our main results are summarized as follows. 

(1) The LMC and SMC become dynamically coupled 
only recently, suffering two strong tidal interactions ~ 2 
Gyr ago and ~ 250 Myr ago. Portions of the SMC disk 
are stripped away during each of these strong encoun- 
ters, with the MS and LA being created during the first 
encounter, and the Bridge during the second. The tim- 
ing of these two tidal interactions is consistent with the 
correlated star formation epochs determined by Harris & 
Zaritsky (2009). 

(2) The morphology of the MS is well reproduced, in- 
cluding its bifurcation into dual filaments (Putman et 
al. 2003a) and its long extent (Nidever et al. 2010). 
The predicted on-sky location of the filaments provides 
a very strong match to observations, as does their cross- 
ing point. The kinematics of the MS are also well repro- 
duced, particularly the velocity inflection at the MS tip 
(Nidever et al. 2010). 

(3) We interpret the observed MS bifurcation as a 
wrapping of the original tidal arm ( "tidal arm A" ) that 
was stripped from the SMC disk. The dynamical process 
responsible for this "wrapping" can be traced to a com- 
plex three-body interplay between the SMC, LMC, and 
MW over the past ~ 2 Gyr. We accordingly suggest that 
the HI morphology of the MS places strong constraints 
on the orbital history of the MCs, namely that the MCs 
have executed at least two pericentric passages about the 
MW during a ~ 2 Gyr bound association. 

(4) The formation of the Bridge ~ 250 Myr ago was 
accompanied by the formation of a complementary struc- 
ture called the "Counter-Bridge" that extends ~ 20 kpc 
along the line of sight behind the SMC. Because of this 
alignment, the Counter-Bridge does not appear as an 



identifiable structure in column density maps, but it is 
likely traced by the large line-of-sight depth of the SMC 
(e.g., Crowl et al. 2001; Groenewegen 2000). 

(5) Two branches of leading material are stripped from 
the SMC disk, but the on-sky location and kinematics do 
not correspond well with the observed properties of the 
LA. Even though we do not explicitly account for hy- 
drodynamical effects in our simulation, we have argued 
that the evolution of the leading material over the past 
~ 300 Myr should be strongly influenced by rampres- 
sure from the hot halo of the Milky Way (Section 5.2 and 



Fig 19 1. Accordingly, we suggest that gravity alone can 



not explain the present location and kinematics of the 
LA, even if this structure is tidal in origin. 

(6) The LMC engulfs a large fraction of the material 
that is tidally stripped from the SMC following their 
strong tidal interactions ~ 2 Gyr ago and ~ 250 Myr 
ago. The transfer of SMC disk material to the LMC 
would correspond to the infall of metal-poor gas and 
would be imprinted into the chemical enrichment his- 
tory of the LMC (e.g., Fig 1 of van Loon et al. 2005; 
Fig 20 of Harris & Zaritsky 2009; Rubele et al. 2012) as 
well as the nitrogen-poor stellar populations of the LMC 
(e.g., Bekki & Tsujimoto 2010). The mass transfer from 
the SMC spheroid would correspond to the accretion of 
old-type stars onto the LMC and is suggestive of vari- 
ous observed populations of kinematically peculiar stars 
(e.g., Kunkel et al. 1997a, Olsen et al. 2011). 

(7) Of our three models for the SMC spheroid, we find 
that the extended spheroid is best able to reproduce the 
relevant observations, including the stellar kinematics of 
the SMC (Harris & Zaritsky 2006), the recent discov- 
ery of a break population of red giants (Nidever et al. 
2011), and the observation of kinematically and chemi- 
cally peculiar stars within the LMC (Olsen et al. 2011). 
Considering these results, future observations may need 
to readdress the possibility that SMC stars were tidally 
stripped into the Bridge region. 

(8) The tidal evolution of the extended spheroid pre- 
dicts that a stellar stream should be stripped away ~ 2 
Gyr ago at the same time that the MS is stripped 
from the SMC disk. Whereas the gaseous MS origi- 
nates in the rotating SMC disk, the stellar stream origi- 
nates in the pressure-supported spheroidal component of 
the SMC. This distinct origin suggests that the stripped 
stars should not correlate with the peak HI regions of 
the MS, as observed (e.g., Recillas-Cruz 1982; Bruck & 
Hawkins 1983; Guhathakurta & Reitzel 1998). Future 
wide-field surveys such as the Southern Sky Survey by 
the SkyMapper telescope (Keller et al. 2007) would be 
ideal for detecting this possible tidal stream and estab- 
lishing its connection to the stellar halo of the SMC. 
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APPENDIX 



In this paper we have presented an N-body model for the Magellanic system characterized by the set of parameter 
values given in Tables 1 and 2. Here we briefly discuss parameter choices other than those adopted to give a sense 
of the parameter space that was explored. In particular, we will discuss Galactic parameters, proper motions of the 
MCs, masses of the MCs, disk angles for the SMC, and the initial time of the simulation. 

Panel (a) of Figure [2ffl giv es the on-sky projection of the "alternate model" as described in Sec 2.1 realized as a 
pure-disk N-body model (Sec |2.3.1 ). The parameter values of this N-body model are the same as those of the adopted 
model except for a different characterization of the Milky Way's NFW halo (see Fig [2] and footnote of Table [IJ). These 
two models share a number of properties, including the fact that their tidally stripped material has nea rly identical 
morphology and kinematics. This similarity in morphology can be assessed by comparing panel (a) of Fig [20^ with the 
left-hand panel of Fig [5] Of particular importance is that the MS exhibits a bifurcation into two distinct filaments 
in both models. The overall kinematics of the MS including the velocity inflection (Nidever et al. 2010) are also 
reproduced in the alternate model, though we do not show the kinematics here. 

The alternate model highlights a strength of the present study, namely that the formation of the MS is consistent 
with multiple prescriptions for the Milky Way halo. The model is therefore more robust in the face of observational 
uncertainties for a number of Galactic parameters, including the mass of the Milky Way (e.g., Gnedin et al. 2010) and 
the circular velocity V c i T at the solar radius (e.g., Kerr & Lynden-Bell 1986; Reid et al. 2009). For instance, we find 
that values of V C i r —240 km/s (adopted model) and V c i r =220 km/s (alternate model) are consistent with the formation 
of the MS without prejudice or preference for either value. Though we found only two possible models in our work, 
it is likely that there exist many other parameterizations of the Milky Way's halo that can generate similar models 
for the MS. Furthermore, it would be worthwhile to consider whether other halo profiles (i.e., other than NFW) can 
generate similar models. 

The present model is less flexible, however, when considering different pr ope r motion values for the LMC and SMC, 
especially those of Kalliva yali l et al. (2006a,b; K06). Panel (b) of Figure 20 provides the on-sky distribution of two 
test particle models (Sec |2.2| _ using the proper motions of Diaz & Bekki(2011a), which is within ler of the K06 
measured values (see Figure [T]) . The black (turquoise) particles are evolved to the present day from an undisturbed 
SMC disk at t = —2.5 Gyr (t = —5.0 Gyr), and the orbital separations are given as the solid (dashed and solid) 
lines in Figure [2l] As is clear, both cases (black and turquoise) give poor results. The stripped material does not 
provide a morphological match to any of the HI structures of the Magellanic system, and the kinematics (not shown) 
are similarly unacceptable. The inability of our adopted model to accommodate the proper motions of K 06 is perhaps 
unsurprising, since the orbital interaction history changes dramatically (e.g., compare Fig [3] and Fig [2]} owing to the 
substitution of much larger proper motions. 

Our decision to establish our parameter space around the measured proper motion values of Vieira et al. (2010, 
V10) rather than those of K06 was informed by a number of factors, including the failure of K06 test particle models 
similar to those shown in panel (b) of Fig 20 Moreover, as mentioned in Sec [II the V10 study possess a number of 
advantages over K06, including a larger sample size, longer baseline, and the ability to measure relative proper motions 
between the LMC and SMC. Lastly, our preference for V10 was also led by the overall properties of the orbital models. 



For instance, Figure 22 gives the fraction of orbits (within the full set of 8 x 10 ° ex plored) in which the SMC suffers 



exactly two close encounters with the LMC. As emphasized below (see also Sec 2.1), there is mounting evidence that 
such a condition on the orbital history is preferred for the formation of the MS (e.g., Besla et al. 2010; Ruzicka et 
al. 2010; Dia z fc Bekki 2 011a ; present work). Of course, there are many more constraints that must be considered 
(e.g., see Sec |2.1| and Sec 2.2), but the condition of exactly two LMC-SMC encounters is a good starting point for 
considering which regions ot parameter space are promising for MS formation scenarios. Fig [22] clearly shows that the 
proper motions with the largest fraction of promising orbits are much closer to the V10 valuesthan those of K06. This 
implies that a parameter space centered on V10 may be more likely to harbor a successful formation scenario for the 
MS, as in the present work. This conclusion is bolstered by the work of Ruzicka et al. (2010), who could not find a 
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f o =-5.0 Gyr 
• fo—2.5 Gyr 
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Fig. 2 0. — The disruption of the SMC disk with different parameter choices, (a) On-sky projection of the "alternate model" (see Tablc[T] 
and Sec |2.1|l w hich differs from the adopted model in the choice of parameters for the Milky Way halo. The N-body model is a pure 
disk (Sec |2.3Tl[ l h aving N = 10 5 total particles and disk angles 9^ = —40° and (f)^ = 230°. (b) On-sky projection of test particle models 
(N = f0 5 ; see Sec \2.2\ with proper motions taken from Diaz & Bekki (20ffa), which is consistent with the measurements of Kallivayalil 
et al. (2006a, b). All other parameters are taken from Table 111 except for the mass of the LMC, which is taken to be 2 X 10 10 M as in 
Diaz & Bekki (2011a). The test particle models differ in theTlme at which each simulation is begun, t = —2.5 Gyr (black) or t = —5.0 
Gyr (turquoise). Coordinates in both panels are galactic longitude I (straight lines) and latitude b (concentric circles) shown in a ZEA 
projection centered on the South Galactic Pole. The orbit trails of the SMC (green) and LMC (red) are also shown, as are their current 
locations (circles). 
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Fig. 21.— Orbital separations between the LMC and MW (red), SMC and MW (green), and LMC and SMC (black) for the two test 
particle models given in panel (b) of Figure [20| Solid lines indicate the time window from t = —2.5 Gyr to the present day, whereas the 
dotted lines begin at an initial time of t = —5.0 Gyr. 



compelling model for the MS within a parameter space centered on the K06 proper motions. I t is also striking to note 
that the proper motions of our adopted model (Fig [T]) are close to the peak regions in Figure 22 



Despite the conclusions of the present work, it may still be possible to reconcile the KOo values with a viable 
MS formation scenario. Furthermore, even though we could not find a first passage model that could explain the 
observations of the MS, the first passage scenario is not necessarily ruled out. It would seem that a different set of 
assumptions (i.e., differing from those of Sec [2]) would be needed to build a model with the K06 proper motions. 
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Fig. 22. — The fraction of orbits having exactly two encounters between the LMC and SMC, plotted as a function of the proper motions 
of the LMC (left) and SMC (right). In each panel, we consider all ~ 8 X 10 6 orbits computed in the present study, and the fraction is 
relative to the total number of orbits computed at each proper motion value. For example, at each point in the LMC proper motion plane, 
there are ~ 5 X 10 4 total orbital models, and the color indicates the fraction of these which satisfy the condition of two encounters. An 
"encounter" is defined as an approach of less than 10 kpc that is also a minimum in the LMC-SMC separation curve (e.g., Figpj] Fig |21[ > . 
See discussion in the Appendix indicating that the condition of two LMC-SMC encounters is a promising criterion for the formation of trie 
MS. Observational ellipses outline the 68.3% confidence regions and lc error bars for Vieira et al. (2010; V10, solid) and Kallivayalil et al. 
(2006a,b; K06, dashed). 

As examples, consider the massive isothermal halo for the Milky Way adopted by Diaz & Bekki (2011a) or efficient 
dynamical friction between the MCs utilized by Besla et al. (2010; see also Besla et al. 2012). One would also need to 
investigate different masses for the LMC and SMC. For instance, Besla et al. (2010) assume that the LMC and SMC 
have halo masses of 1.8xlO n M and 2.5xl0 10 M , respectively, which are roughly an order of magnitude larger than 
adopted in the present study (Table [T]) . The reason behind this large difference in adopted MC masses may possibly 
be traced to the difference in model assumptions between the two studies. In particular, mutual dynamical friction 
between the MCs is ignored in the present study, but in Besla et al. (2010), it is the mechanism that drives the MCs 
toward the mutual encounter that creates the MS. It is not surprising that a different range of masses for the LMC 
and SMC are preferred in the presence/absence of dynamical friction, especially because the magnitude of the effect 
is sensitive to satellite mass (e.g., see equation (JsJ) ) . 

In the present work, we explored a range of masses similar to those of traditional tidal models (e.g., GN96, C06): (1.0, 
1.5, 2.0, 2.5, 3.5, 4.0) xlO 10 M Q for the LMC, and (3.0, 4.5) xlO 9 M for the SMC. We could find promising formation 
scenarios for the MS with only one mass combination (i.e., l.OxlO 10 M and 3.0xl0 9 M ), which certainly falls in 
the "low mass" regime for the MCs (e.g., van der Marel et al. 2002, Stanimirovic et al. 2004, Bekki & Stanimirovic 
2009). This may be a reflection of the actual masses of the MCs (or possibly their mass ratio), but one must not forget 
that the adopted values are a product of our particular model assumptions, just as in the case of large MC masses for 
Besla et al. (2010). Among our assumptions is the fact that the LMC is represented by a fixed Plummer potential 



(Sec 2.1 ), which underpins the lack of dynamical friction between the MCs. For the LMC itself, this assumption is not 

■>2 Gyr, as verified by live models of the LMC (e.g. 



unreasonable since its mass is unlikely to change within the past 

Bekki & Chiba 2005). It remains unclear, however, if our neglect of dynamical friction between the MCs is realistic. 

F igure 23 shows how our N-body model is affected by adjusting the SMC disk angles 9d and (f>d (see definitions in 
Sec |2.3.1 )~Tn general, we found that the best disk angles fell in the range -60° < 8 d < -40° and 220° < 4>d < 260°, 
and models outside of this range failed to reproduce a tidal stream resembling the MS. Even within thi s ra nge, however, 

dashed line) 



23 



many models provide a poor match to the MS, with either too many particles at the MS tip (Fig 
or too few particles in the MS altogether (dotted). The reasons for this sensitivity to disk orientation have not yet 
been fully explored, nor is it clear why there seems to be a preference for = —45° and <f>d = 230° regardless of tidal 
history (i.e, these values were found to be best in GN96 as well as the present work). Although the density profile 
of the MS in our adopted model (Fig 23 solid line) is superior to other candidate models, there are discrepancies 
with observation. The column density of the MS is observed to decline exponentially along its length (e.g., Figure 
10 of Nidever et al. 2010), whereas the particle count in the adopted model exhibits a small peak at the MS tip 
Lms ~ —120°, and furthermore the density along the MS does not decli ne a s rapidly as observed. Reconciling theory 
with observation will likely require the inclusion of hydrodynamics (Sec 5.2). 



The final input parameter that we discuss is the initial time of the simulation. Despite the poor observational 
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Fig. 23. — Particle count as a function of Magellanic longitude Ljyis quoted as a percentage of the total for three models: the disk of 
Model 1 (solid 6» d = -45° d = 230°) and two different pure disk models (dashed G d = -40° cf> d = 260°; and dotted 8 d = -60° <^ d = 230°). 
See Sec |2.3.i"1 for the definition of pure disk models. The density profile of the MS in each model is given at Lms < —25°, and that of the 
LA at L M s > 0° . 

the model that evolves from t = —2.5 Gyr (black) 
that evolves from t = —5.0 Gyr (turquoise) whose 



correspondence of the test particle models in panel (b) of Fig 20 
is able to form coherent tidal structures in contrast to the moae 
str uctu res are "smeared" spatially and kinematically. This difference can be traced to the interaction history given in 
Fig [2l] The SMC suffers two close encounters with the LMC in the black model (solid lines), similar to the scenario 
in the present work. The turquoise model evolves through these interactions as well, but only after suffering earlier 
encounters with the LMC and MW between t = —5.0 Gyr and t — —2.5 Gyr (dashed lines). The multiplicity of these 
interactions effectively destroys the coherence of the stripped material. Similarly, GSF94 find that the coherence of 
the MS is destroyed if the SMC is subjected to tidal interactions previous to ~3 Gyr ago. Based on these results and 
those of other recent models (e.g. Besla et al. 2010; Ruzicka et al. 2010; Diaz & Bekki 2011a), we suggest that the 
creation of a coherent tidal structure resembling the observed MS requires exactly two recent encounters between the 
LMC and SMC. 

The time at which the simulation is begun is therefore an important parameter for models having earlier epochs of 
strong tidal interaction (e.g., GSF94, GN96, C06). For example, if the LMC and SMC are indefinitely bound to one 



another (as in the models of panel (b) in Fig 20 ) , one may artificially ensure that only two strong interactions occur 
by appropriately choosing the initial time for the simulation. This tactic is unnecessary for models that assume the 
LMC and SMC have only recently entered their bound state (Besla et al. 2010; Diaz & Bekki 2011a; present work). 
A recent dynamical coupling between the LMC and SMC can elegantly deliver the two necessary tidal encounters for 
the formation of the MS while keeping the SMC safe from external tidal fields at earlier epochs. In principle, the 
starting time of the simulation should not matter for such models as long as it is prior to the first strong encounter. 
We mention briefly a caveat of the present model, however, because the initial time of our adopted model is chosen 
at t = —3.37 Gyr. We found that the internal structure of the MS depends somewhat on the orientation of the SMC 
bar at the first epoch of tidal stripping, and since the phase of the SMC bar depends on the initial start time of the 
simulation, the present results depend weakly on the initial time. 



